Re: Use of C99 extra long double math functions after r236148
Am 26.07.2012 03:52 (UTC+1) schrieb Bruce Evans: On Wed, 25 Jul 2012, Rainer Hurling wrote: On 25.07.2012 19:00 (UTC+2), Steve Kargl wrote: On Wed, Jul 25, 2012 at 06:29:18PM +0200, Rainer Hurling wrote: Many thanks to you three for implementing expl() with r238722 and r238724. I am not a C programmer, but would like to ask if the following example is correct and suituable as a minimalistic test of this new C99 function? It's not clear to me what you mean by test. If expl() is not available in libm, then linking the code would fail. So, testing for the existence of expl() (for example, in a configure script) is as simple as Sorry for not being clear enough. I didn't mean testing for the existence, but for some comparable output between exp() and expl(), on a system with expl() available in libm. This is basically what I do to test exp() (with a few billion cases automatically generated and compared). It is not sufficient for checking expl(), except for consistency. (It is assumed that expl() is reasonably accurate. If it is in fact less accurate than exp(), this tends to show up in the comparisons.) Ahh ok, I think I understand it better now. #include math.h long double func(long double x) { return (expl(x)); } //--- #include stdio.h #include math.h int main(void) { double c = 2.0; long double d = 2.0; double e = exp(c); long double f = expl(d); printf(exp(%f) is %.*f\n, c, 90, e); printf(expl(%Lf) is %.*Lf\n, d, 90, f); If you mean testing that the output is correct, then asking for 90 digits is of little use. The following is sufficient (and my actually produce a digit or two more than is available in number) Ok, I understand. I printed the 90 digits to be able to take a look at the decimal places, I did not expect to get valid digits in this area. Use binary format (%a) for manual comparison. Don't print any more bits than the format has. This is DBL_MANT_DIG (53) for doubles and LDLBL_MANT_DIG (64 on x86) for long doubles. %a format is in nybbles and tends to group the bits into nybbles badly. See below on reducing problems from this. Decimal format has to print about 3 more digits than are really meaningful, to allow recovering the original value uniquely. For manual comparison, you need to print these extra digits and manually round or ignore them as appropriate. The correct number of extra digits is hard to determine. For the any, type, it is DECIMAL_DIG (21) on x86. The corresponding number of normally-accurate decimal digits for long doubles is given by LDBL_DIG (18). For floats and doubles, this corresponds to FLT_DIG (6) and DBL_DIG (15). Unfortunately, float.h doesn't define anything corresponding to DECIMAL_DIG for the smaller types. 21 is a lot of digits and noise digits take a long time to determine and ignore (its worse on sparc64 where DECIMAL_DIG is 36). I usually add 2 extra digits to the number of normally-accurate digits. This is sloppy. 3 is needed in some cases, depending on MANT_DIG and the bits in log(2) and/or log(10). I printed more bits (digits) than the format provides, because I wanted to see if and when what the new function would do in this 'outside' area. Of course you are right that this has nothing to do with more precision and these digits may not be used in subsequent calculations. But this really was not my intention here (only curiosity). troutmask:fvwm:kargl[203] diff -u a.c.orig a.c --- a.c.orig2012-07-25 09:38:31.0 -0700 +++ a.c 2012-07-25 09:40:36.0 -0700 @@ -1,5 +1,6 @@ #include stdio.h #include math.h +#include float.h int main(void) { @@ -9,8 +10,8 @@ double e = exp(c); long double f = expl(d); - printf(exp(%f) is %.*f\n, c, 90, e); - printf(expl(%Lf) is %.*Lf\n, d, 90, f); + printf(exp(%f) is %.*f\n, c, DBL_DIG+2, e); + printf(expl(%Lf) is %.*Lf\n, d, LDBL_DIG+2, f); return 0; } Thanks, I was not aware of DBL_DIG and LDBL_DIG. Steve is sloppy and adds 2 also :-). For long doubles, it is clear that 3 are strictly needed, since DECIMAL_DIG is 3 more. So I better use +3 here? printf(expl(%Lf) is %.*Lf\n, d, LDBL_DIG+3, f); ^^ For most long double functions on i386, you need to switch the rounding precision to 64 bits around calls to them, and also to do any operations on the results except printing them. expl() is one of the few large functions that does the switch internally. So the above should work (since it only prints), but (expl(d) + 0) should round to the default 53-bit precision and this give the same result as exp(d). Now that this rounding problem is more understandable to me, I am happy that in most cases I can use math/R for my calculations ;-) If you actually want to test expl() to see if it is producing a decent result, you need a reference solution that contains a higher precision. I use mpfr with 256 bits of precision.
Re: Use of C99 extra long double math functions after r236148
On Wed, Jul 25, 2012 at 06:29:18PM +0200, Rainer Hurling wrote: Many thanks to you three for implementing expl() with r238722 and r238724. I am not a C programmer, but would like to ask if the following example is correct and suituable as a minimalistic test of this new C99 function? It's not clear to me what you mean by test. If expl() is not available in libm, then linking the code would fail. So, testing for the existence of expl() (for example, in a configure script) is as simple as #include math.h long double func(long double x) { return (expl(x)); } //--- #include stdio.h #include math.h int main(void) { double c = 2.0; long double d = 2.0; double e = exp(c); long double f = expl(d); printf(exp(%f) is %.*f\n, c, 90, e); printf(expl(%Lf) is %.*Lf\n, d, 90, f); If you mean testing that the output is correct, then asking for 90 digits is of little use. The following is sufficient (and my actually produce a digit or two more than is available in number) troutmask:fvwm:kargl[203] diff -u a.c.orig a.c --- a.c.orig2012-07-25 09:38:31.0 -0700 +++ a.c 2012-07-25 09:40:36.0 -0700 @@ -1,5 +1,6 @@ #include stdio.h #include math.h +#include float.h int main(void) { @@ -9,8 +10,8 @@ double e = exp(c); long double f = expl(d); - printf(exp(%f) is %.*f\n, c, 90, e); - printf(expl(%Lf) is %.*Lf\n, d, 90, f); + printf(exp(%f) is %.*f\n, c, DBL_DIG+2, e); + printf(expl(%Lf) is %.*Lf\n, d, LDBL_DIG+2, f); return 0; } return 0; } //--- Compiled with 'c99 -o math_expl math_expl.c -lm' and running afterwards it gives me: exp(2.00) is 7.38905609893065040 expl(2.00) is 7.38905609893065022739 Already the sixteenth position after decimal point differs. DBL_DIG is 15 on x86 hardware and LDBL_DIG is 18. You can expect at most 15 correct digits from exp(). Please forgive me, if this is a really stupid approach for producing some expl() output. If you actually want to test expl() to see if it is producing a decent result, you need a reference solution that contains a higher precision. I use mpfr with 256 bits of precision. troutmask:fvwm:kargl[213] ./testl -V 2 ULP = 0.3863 x = 2.00e+00 libm: 7.389056098930650227e+00 0x1.d8e64b8d4ddadcc4p+2 mpfr: 7.389056098930650227e+00 0x1.d8e64b8d4ddadcc4p+2 mpfr: 7.3890560989306502272304274605750078131803155705518\ 47324087127822522573796079054e+00 mpfr: 0x7.63992e35376b730ce8ee881ada2aeea11eb9ebd93c887eb59ed77977d109f148p+0 The 1st 'mpfr:' line is produced after converting the results fof mpfr_exp() to long double. The 2nd 'mpfr:' line is produced by mpfr_printf() where the number of printed digits depends on the 256-bit precision. The last 'mpfr:' line is mpfr_printf()'s hex formatting. Unfortunately, it does not normalize the hex representation to start with '0x1.', which makes comparison somewhat difficult. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 11.07.2012 00:58 (UTC+2), David Schultz wrote: On Tue, Jul 10, 2012, Rainer Hurling wrote: On 10.07.2012 17:11 (UTC+2), David Schultz wrote: On Tue, Jul 10, 2012, Rainer Hurling wrote: On 10.07.2012 16:02 (UTC+2), Warner Losh wrote: On Jul 10, 2012, at 3:10 AM, Rainer Hurling wrote: As far as I understand from discussions on R mailing list (r-de...@r-project.org), they plan to reduce the emulation and/or workaround of long and complex math functions for FreeBSD and other systems with their next releases of R devel. So we could really need some progress with our C99 conform math functions ;-) Not having R would be a bit pain in my backside. That's one of the practical considerations that I was talking about. It is very real, and if I have to, I'll commit the #define junk I railed against to get it back. Please, let's get some progress. I have some time to help. Yes, thank you Warner, that is also my problem. As I wrote some weeks ago (05/28/2012) when starting this thread, I am using FreeBSD as a scientific desktop because of its good scaling properties. For some years now, FreeBSD fits all our needs with R, SAGA GIS, PostgreSQL and some more. If I would not be able to run upcoming versions of R on FreeBSD any more, that would be really, really hard :-( Do you have a list of the essential functions here? There are 17 long double functions and some complex functions missing, but only a handful of those are of general interest. The reason I ask is that if R is just looking for a few missing functions that are already mostly implemented, then the best solution is probably to finish that work. But if it's expecting us to have something arcane like long double Bessel functions of the first kind, then we need to pursue a workaround in the short term. That is, what I found by grepping the sources of a recent R development version: expl: src/nmath/pnchisq.c Many thanks to you three for implementing expl() with r238722 and r238724. I am not a C programmer, but would like to ask if the following example is correct and suituable as a minimalistic test of this new C99 function? //--- #include stdio.h #include math.h int main(void) { double c = 2.0; long double d = 2.0; double e = exp(c); long double f = expl(d); printf(exp(%f) is %.*f\n, c, 90, e); printf(expl(%Lf) is %.*Lf\n, d, 90, f); return 0; } //--- Compiled with 'c99 -o math_expl math_expl.c -lm' and running afterwards it gives me: exp(2.00) is 7.38905609893065040694182243896648287773132324218750 expl(2.00) is 7.38905609893065022739794267536694860609713941812515258789062500 Already the sixteenth position after decimal point differs. Please forgive me, if this is a really stupid approach for producing some expl() output. uname -a: FreeBSD xxx.xxx.xxx 10.0-CURRENT FreeBSD 10.0-CURRENT #0 r238740: Tue Jul 24 18:08:13 CEST 2012 x...@xxx.xxx.xxx:/usr/obj/usr/src/sys/XXX amd64 Rainer logl: src/nmath/dnbeta.c src/nmath/pnbeta.c Bruce has versions of these that could be committed with some cleanup. It's a matter of sorting through about 1200 emails from him and 3 source trees to find the most up-to-date patches, then cleaning them up and testing and committing them. I have no time right now, but I will do at least the first step as soon as I can, and try to get the patches to someone willing to do the final few steps. log10l: src/extra/trio/trio.c log1pl: src/nmath/pnbeta.c If Bruce doesn't already have implementations of these, they are easy wrappers around logl() or some internal k_logl in Bruce's implementation. powl: src/extra/trio/triostr.c src/extra/trio/trio.c src/main/format.c It's hard to do a good job on powl(), but the simple approach (exp(log(x)*y)) plus a few special cases may suffice for many uses. NEWS:l2044 The C99 functions acosh, asinh, atanh, snprintf and vsnprintf are now required. We have had them forever. NEWS:l3032 The C99 double complex type is now required. The C99 complex trigonometric functions (such as csin) are not currently required (FreeBSD lacks most of them): substitutes are used if they are missing. We have these (but not the inverse trig functions). NEWS:l3277 Complex arithmetic (notably z^n for complex z and integer n) gave incorrect results since R 2.10.0 on platforms without C99 complex support. This and some lesser issues in trigonometric functions have been corrected. Such platforms were rare (we know of Cygwin and FreeBSD). However, because of new compiler optimizations in the way complex arguments are handled, the same code was selected on x86_64 Linux with gcc 4.5.x at the default -O2 optimization (but not at -O). Not sure if this is relevant. BTW: There seems to be a discrepancy about missing functions listed in http://wiki.freebsd.org/MissingMathStuff
Re: Use of C99 extra long double math functions after r236148
On 07/25/12 11:29, Rainer Hurling wrote: Many thanks to you three for implementing expl() with r238722 and r238724. I am not a C programmer, but would like to ask if the following example is correct and suituable as a minimalistic test of this new C99 function? //--- #include stdio.h #include math.h int main(void) { double c = 2.0; long double d = 2.0; double e = exp(c); long double f = expl(d); printf(exp(%f) is %.*f\n, c, 90, e); printf(expl(%Lf) is %.*Lf\n, d, 90, f); return 0; } //--- Compiled with 'c99 -o math_expl math_expl.c -lm' and running afterwards it gives me: exp(2.00) is 7.38905609893065040694182243896648287773132324218750 expl(2.00) is 7.38905609893065022739794267536694860609713941812515258789062500 Just as a point of comparison, here is the answer computed using Mathematica: N[Exp[2], 50] 7.3890560989306502272304274605750078131803155705518 As you can see, the expl solution has only a few digits more accuracy that exp. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Wed, Jul 25, 2012 at 12:27:43PM -0500, Stephen Montgomery-Smith wrote: On 07/25/12 11:29, Rainer Hurling wrote: Many thanks to you three for implementing expl() with r238722 and r238724. I am not a C programmer, but would like to ask if the following example is correct and suituable as a minimalistic test of this new C99 function? (program deleted) Compiled with 'c99 -o math_expl math_expl.c -lm' and running afterwards it gives me: exp(2.00) is 7.3890560989306504069 expl(2.00) is 7.38905609893065022739794 Just as a point of comparison, here is the answer computed using Mathematica: N[Exp[2], 50] 7.3890560989306502272304274605750078131803155705518 As you can see, the expl solution has only a few digits more accuracy that exp. Unless you are using sparc64 hardware. flame:kargl[204] ./testl -V 2 ULP = 0.2670 for x = 2.0e+00 mpfr exp: 7.389056098930650227230427460575008e+00 libm exp: 7.389056098930650227230427460575008e+00 -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 25.07.2012 19:00 (UTC+2), Steve Kargl wrote: On Wed, Jul 25, 2012 at 06:29:18PM +0200, Rainer Hurling wrote: Many thanks to you three for implementing expl() with r238722 and r238724. I am not a C programmer, but would like to ask if the following example is correct and suituable as a minimalistic test of this new C99 function? It's not clear to me what you mean by test. If expl() is not available in libm, then linking the code would fail. So, testing for the existence of expl() (for example, in a configure script) is as simple as Sorry for not being clear enough. I didn't mean testing for the existence, but for some comparable output between exp() and expl(), on a system with expl() available in libm. #include math.h long double func(long double x) { return (expl(x)); } //--- #include stdio.h #include math.h int main(void) { double c = 2.0; long double d = 2.0; double e = exp(c); long double f = expl(d); printf(exp(%f) is %.*f\n, c, 90, e); printf(expl(%Lf) is %.*Lf\n, d, 90, f); If you mean testing that the output is correct, then asking for 90 digits is of little use. The following is sufficient (and my actually produce a digit or two more than is available in number) Ok, I understand. I printed the 90 digits to be able to take a look at the decimal places, I did not expect to get valid digits in this area. troutmask:fvwm:kargl[203] diff -u a.c.orig a.c --- a.c.orig2012-07-25 09:38:31.0 -0700 +++ a.c 2012-07-25 09:40:36.0 -0700 @@ -1,5 +1,6 @@ #include stdio.h #include math.h +#include float.h int main(void) { @@ -9,8 +10,8 @@ double e = exp(c); long double f = expl(d); - printf(exp(%f) is %.*f\n, c, 90, e); - printf(expl(%Lf) is %.*Lf\n, d, 90, f); + printf(exp(%f) is %.*f\n, c, DBL_DIG+2, e); + printf(expl(%Lf) is %.*Lf\n, d, LDBL_DIG+2, f); return 0; } Thanks, I was not aware of DBL_DIG and LDBL_DIG. return 0; } //--- Compiled with 'c99 -o math_expl math_expl.c -lm' and running afterwards it gives me: exp(2.00) is 7.38905609893065040 expl(2.00) is 7.38905609893065022739 Already the sixteenth position after decimal point differs. DBL_DIG is 15 on x86 hardware and LDBL_DIG is 18. You can expect at most 15 correct digits from exp(). Please forgive me, if this is a really stupid approach for producing some expl() output. If you actually want to test expl() to see if it is producing a decent result, you need a reference solution that contains a higher precision. I use mpfr with 256 bits of precision. troutmask:fvwm:kargl[213] ./testl -V 2 ULP = 0.3863 x = 2.00e+00 libm: 7.389056098930650227e+00 0x1.d8e64b8d4ddadcc4p+2 mpfr: 7.389056098930650227e+00 0x1.d8e64b8d4ddadcc4p+2 mpfr: 7.3890560989306502272304274605750078131803155705518\ 47324087127822522573796079054e+00 mpfr: 0x7.63992e35376b730ce8ee881ada2aeea11eb9ebd93c887eb59ed77977d109f148p+0 The 1st 'mpfr:' line is produced after converting the results fof mpfr_exp() to long double. The 2nd 'mpfr:' line is produced by mpfr_printf() where the number of printed digits depends on the 256-bit precision. The last 'mpfr:' line is mpfr_printf()'s hex formatting. Unfortunately, it does not normalize the hex representation to start with '0x1.', which makes comparison somewhat difficult. Many thanks also for this mpfr example. This helps me to understand a little bit more what is going here. So on amd64 even the expl() result is far from what mpfr provides. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 07/25/12 12:31, Steve Kargl wrote: On Wed, Jul 25, 2012 at 12:27:43PM -0500, Stephen Montgomery-Smith wrote: On 07/25/12 11:29, Rainer Hurling wrote: Many thanks to you three for implementing expl() with r238722 and r238724. I am not a C programmer, but would like to ask if the following example is correct and suituable as a minimalistic test of this new C99 function? (program deleted) Compiled with 'c99 -o math_expl math_expl.c -lm' and running afterwards it gives me: exp(2.00) is 7.3890560989306504069 expl(2.00) is 7.38905609893065022739794 Just as a point of comparison, here is the answer computed using Mathematica: N[Exp[2], 50] 7.3890560989306502272304274605750078131803155705518 As you can see, the expl solution has only a few digits more accuracy that exp. Unless you are using sparc64 hardware. flame:kargl[204] ./testl -V 2 ULP = 0.2670 for x = 2.0e+00 mpfr exp: 7.389056098930650227230427460575008e+00 libm exp: 7.389056098930650227230427460575008e+00 Yes. It would be nice if long on the Intel was as long as the sparc64. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Wed, Jul 25, 2012 at 07:33:07PM +0200, Rainer Hurling wrote: On 25.07.2012 19:00 (UTC+2), Steve Kargl wrote: If you actually want to test expl() to see if it is producing a decent result, you need a reference solution that contains a higher precision. I use mpfr with 256 bits of precision. troutmask:fvwm:kargl[213] ./testl -V 2 ULP = 0.3863 x = 2.00e+00 libm: 7.389056098930650227e+00 0x1.d8e64b8d4ddadcc4p+2 mpfr: 7.389056098930650227e+00 0x1.d8e64b8d4ddadcc4p+2 mpfr: 7.3890560989306502272304274605750078131803155705518\ 47324087127822522573796079054e+00 mpfr: 0x7.63992e35376b730ce8ee881ada2aeea11eb9ebd93c887eb59ed77977d109f148p+0 The 1st 'mpfr:' line is produced after converting the results fof mpfr_exp() to long double. The 2nd 'mpfr:' line is produced by mpfr_printf() where the number of printed digits depends on the 256-bit precision. The last 'mpfr:' line is mpfr_printf()'s hex formatting. Unfortunately, it does not normalize the hex representation to start with '0x1.', which makes comparison somewhat difficult. Many thanks also for this mpfr example. This helps me to understand a little bit more what is going here. So on amd64 even the expl() result is far from what mpfr provides. Of course!. MPFR is a multiple precision library. One specifies the precision, and mpfr returns a value with that precision. #include mpfr.h int main(void) { int i, j[5] = {24, 53, 64, 113, 256}; mpfr_t x, f; for (i = 0; i 5; i++) { /* Set working precision to j[i]. */ mpfr_inits2(j[i], x, f, MPFR_RNDN); mpfr_set_ui(x, 2, MPFR_RNDN); mpfr_exp(f, x, MPFR_RNDN); mpfr_printf(exp(%Re) = %Re\n, x, f); mpfr_clears(x, f, NULL); } } troutmask:fvwm:kargl[222] cc -o z -I/usr/local/include a.c -L/usr/local/lib\ -lmpfr -lgmp troutmask:fvwm:kargl[223] ./z exp(2e+00) = 7.38905621e+00 exp(2e+00) = 7.3890560989306504e+00 exp(2e+00) = 7.3890560989306502274e+00 exp(2e+00) = 7.38905609893065022723042746057500802e+00 exp(2e+00) = 7.38905609893065022723042746057500781\ 3180315570551847324087127822522573796079054e+00 -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Wed, 25 Jul 2012, Rainer Hurling wrote: On 25.07.2012 19:00 (UTC+2), Steve Kargl wrote: On Wed, Jul 25, 2012 at 06:29:18PM +0200, Rainer Hurling wrote: Many thanks to you three for implementing expl() with r238722 and r238724. I am not a C programmer, but would like to ask if the following example is correct and suituable as a minimalistic test of this new C99 function? It's not clear to me what you mean by test. If expl() is not available in libm, then linking the code would fail. So, testing for the existence of expl() (for example, in a configure script) is as simple as Sorry for not being clear enough. I didn't mean testing for the existence, but for some comparable output between exp() and expl(), on a system with expl() available in libm. This is basically what I do to test exp() (with a few billion cases automatically generated and compared). It is not sufficient for checking expl(), except for consistency. (It is assumed that expl() is reasonably accurate. If it is in fact less accurate than exp(), this tends to show up in the comparisons.) #include math.h long double func(long double x) { return (expl(x)); } //--- #include stdio.h #include math.h int main(void) { double c = 2.0; long double d = 2.0; double e = exp(c); long double f = expl(d); printf(exp(%f) is %.*f\n, c, 90, e); printf(expl(%Lf) is %.*Lf\n, d, 90, f); If you mean testing that the output is correct, then asking for 90 digits is of little use. The following is sufficient (and my actually produce a digit or two more than is available in number) Ok, I understand. I printed the 90 digits to be able to take a look at the decimal places, I did not expect to get valid digits in this area. Use binary format (%a) for manual comparison. Don't print any more bits than the format has. This is DBL_MANT_DIG (53) for doubles and LDLBL_MANT_DIG (64 on x86) for long doubles. %a format is in nybbles and tends to group the bits into nybbles badly. See below on reducing problems from this. Decimal format has to print about 3 more digits than are really meaningful, to allow recovering the original value uniquely. For manual comparison, you need to print these extra digits and manually round or ignore them as appropriate. The correct number of extra digits is hard to determine. For the any, type, it is DECIMAL_DIG (21) on x86. The corresponding number of normally-accurate decimal digits for long doubles is given by LDBL_DIG (18). For floats and doubles, this corresponds to FLT_DIG (6) and DBL_DIG (15). Unfortunately, float.h doesn't define anything corresponding to DECIMAL_DIG for the smaller types. 21 is a lot of digits and noise digits take a long time to determine and ignore (its worse on sparc64 where DECIMAL_DIG is 36). I usually add 2 extra digits to the number of normally-accurate digits. This is sloppy. 3 is needed in some cases, depending on MANT_DIG and the bits in log(2) and/or log(10). troutmask:fvwm:kargl[203] diff -u a.c.orig a.c --- a.c.orig2012-07-25 09:38:31.0 -0700 +++ a.c 2012-07-25 09:40:36.0 -0700 @@ -1,5 +1,6 @@ #include stdio.h #include math.h +#include float.h int main(void) { @@ -9,8 +10,8 @@ double e = exp(c); long double f = expl(d); - printf(exp(%f) is %.*f\n, c, 90, e); - printf(expl(%Lf) is %.*Lf\n, d, 90, f); + printf(exp(%f) is %.*f\n, c, DBL_DIG+2, e); + printf(expl(%Lf) is %.*Lf\n, d, LDBL_DIG+2, f); return 0; } Thanks, I was not aware of DBL_DIG and LDBL_DIG. Steve is sloppy and adds 2 also :-). For long doubles, it is clear that 3 are strictly needed, since DECIMAL_DIG is 3 more. For most long double functions on i386, you need to switch the rounding precision to 64 bits around calls to them, and also to do any operations on the results except printing them. expl() is one of the few large functions that does the switch internally. So the above should work (since it only prints), but (expl(d) + 0) should round to the default 53-bit precision and this give the same result as exp(d). If you actually want to test expl() to see if it is producing a decent result, you need a reference solution that contains a higher precision. I use mpfr with 256 bits of precision. troutmask:fvwm:kargl[213] ./testl -V 2 ULP = 0.3863 x = 2.00e+00 libm: 7.389056098930650227e+00 0x1.d8e64b8d4ddadcc4p+2 mpfr: 7.389056098930650227e+00 0x1.d8e64b8d4ddadcc4p+2 mpfr: 7.3890560989306502272304274605750078131803155705518\ 47324087127822522573796079054e+00 mpfr: 0x7.63992e35376b730ce8ee881ada2aeea11eb9ebd93c887eb59ed77977d109f148p+0 The 1st 'mpfr:' line is produced after converting the results fof mpfr_exp() to long double. The 2nd 'mpfr:' line is produced by mpfr_printf() where the number of printed digits depends on the 256-bit precision. The last 'mpfr:' line is mpfr_printf()'s hex formatting. Unfortunately, it does not normalize the hex
Re: Use of C99 extra long double math functions after r236148
On Wed, 25 Jul 2012, Stephen Montgomery-Smith wrote: On 07/25/12 12:31, Steve Kargl wrote: On Wed, Jul 25, 2012 at 12:27:43PM -0500, Stephen Montgomery-Smith wrote: Just as a point of comparison, here is the answer computed using Mathematica: N[Exp[2], 50] 7.3890560989306502272304274605750078131803155705518 As you can see, the expl solution has only a few digits more accuracy that exp. Unless you are using sparc64 hardware. flame:kargl[204] ./testl -V 2 ULP = 0.2670 for x = 2.0e+00 mpfr exp: 7.389056098930650227230427460575008e+00 libm exp: 7.389056098930650227230427460575008e+00 Yes. It would be nice if long on the Intel was as long as the sparc64. You want it to be as slow as sparc64? (About 300 times slower, after scaling the CPU clock rates. Doubles on sparc64 are less than 2 times slower.) I forgot to mention in a previous reply is that expl has only a few more decimal digits of accuracy than exp because the extra precision on x86 wasn't designed to give much more accuracy. It was designed to give more chance of full double precision accuracy in naive code. It was designed in ~1980 when bits were expensive and the extra 11 provided by the 8087 were considered the best tradeoff between cost and accuracy. They only previde 2-3 extra decimal digits of accuracy. They are best thought of as guard bits. Floating point uses 1 or 2 guard bits internally. 11 extends that significantly and externalizes it, but is far from doubling the number of bits. Their use to provide extra precision was mostly defeated in C by bad C bindings and implementations. This was consolidated by my not using the extra bits for the default rounding precision in FreeBSD. This has been further consolidated by SSE not supporting extended precision. Now the naive code that uses doubles never gets the extra precision on amd64. Mixing of long doubles with doubles is much slower with SSE+i387 than with i387, since the long doubles are handled in different registers and must be translated with SSE+i387, while with i387, using long doubles is almost free (it actually has a negative cost in non-naive code since it allows avoiding extra precision in software). Thus SSE also inhibits using the extra precision intentionally. Bruce ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Fri, Jul 13, 2012, Warner Losh wrote: Just to jump back into the fray a bit, since this point hasn't been articulated well. On Jul 10, 2012, at 6:55 PM, Peter Jeremy wrote: On 2012-Jul-08 19:01:07 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: Well, on the most popular hardware (that being i386/amd64), ld80 will use hardware fp instruction while ld128 must be done completely in software. The speed difference is significant. AFAIK, of the architectures that FreeBSD supports, only sparc64 defines ld128 in the architecture and I don't believe there are any SPARC chip implementations that implement ld128 math in hardware. We shouldn't be gating the new math on an issue that only affects sparc64 machines. If they have ld80 level of support for that architecture, then that is sufficient to get things into the tree. There's no real benefit from making numerics good on sparc64 for the project, since our support for the platform isn't stellar and the platform itself is getting a bit long in the tooth. That said, if people want to do it, be my guest. If it is important enough to catch someone's attention, then it is important enough to have. It just isn't important enough to be a gating factor if nobody has signed up for it yet. I have generally encouraged people to develop both at the same time, for three reasons: 1. Development is more efficient that way. When the algorithm is fresh in your mind, it's just a question of changing the constants and polynomials. 2. The 128-bit format is increasingly being supported on other platforms via software emulation (e.g., __float128 in gcc) and may be more widely supported in hardware in the future. 3. If the ld128 implementations don't happen at the same time as the ld80 versions, it becomes a pain for ports folks, and furthermore, it may *never* get implemented. But you are right that ld128 should not be a blocking issue. If ld128 is preventing people from making progress, then that work should definitely be deferred. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 13 Jul 2012, at 17:40, Warner Losh wrote: We shouldn't be gating the new math on an issue that only affects sparc64 machines Mostly agreed, but it's worth noting that the APCS for ARMv8[1] specifies that long double should be IEEE 754- 2008 quad precision. I think ARMv8 is going to be an important platform for us in the next few years, so ld128 is not just for SPARC (I'd also like us to switch to ld128 on PowerPC, since that's the ABI that everyone else uses, but I don't anticipate PowerPC becoming tier 1 any time soon, whereas I would very much like to see ARMv8 become tier 1 within a year of shipping silicon). That's not to say that we should hold things up waiting for ld128 versions, just that adding ld128 versions soon after adding ld80 ones would be very helpful. David [1] http://infocenter.arm.com/help/topic/com.arm.doc.ihi0055a/IHI0055A_aapcs64.pdf___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Jul 14, 2012, at 3:24 AM, David Chisnall wrote: On 13 Jul 2012, at 17:40, Warner Losh wrote: We shouldn't be gating the new math on an issue that only affects sparc64 machines Mostly agreed, but it's worth noting that the APCS for ARMv8[1] specifies that long double should be IEEE 754- 2008 quad precision. I think ARMv8 is going to be an important platform for us in the next few years, so ld128 is not just for SPARC (I'd also like us to switch to ld128 on PowerPC, since that's the ABI that everyone else uses, but I don't anticipate PowerPC becoming tier 1 any time soon, whereas I would very much like to see ARMv8 become tier 1 within a year of shipping silicon). That's not to say that we should hold things up waiting for ld128 versions, just that adding ld128 versions soon after adding ld80 ones would be very helpful. Agreed. My point was more don't gate everything on the 128-bit versions not hey, nobody do this, it is stupid. If the 128-bit versions are ready at the same time as the 80-bit versions, so much the better. Warner David [1] http://infocenter.arm.com/help/topic/com.arm.doc.ihi0055a/IHI0055A_aapcs64.pdf ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 2012-Jul-11 15:32:47 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: I know an approach to implementing many of the missing functions. Are you willing to share this insight so someone else could do the work? When I do find some free time, I look at what is missing and start to put together a new function. At the moment, it seems that it takes 3+ years to get a new function written, tested, and committed. And, from what I can see, much of this is done quietly - which opens up the possibility that two people might both implement the same code or that people will avoid the area in fear of treading on someone else's toes. As I said previously, I believe the existing wiki page could be improved to form a central co-ordinating point to show what what activity is (or isn't) occurring. but most people seem to push the easy button and want to grab either cephes or netlib's libm. There are technical issues with this approach that I won't rehash again. Doing it properly requires significant effort by people with fairly specialised skills. Whilst the project has several people with the skills, it appears that none of them currently have the time. In the meantime, FreeBSD is taking free kicks from other FOSS groups that have gone down the quick-and-dirty path. AFAIK, none of the relevant standards (POSIX, IEEE754) have any precision requirements for functions other than +-*/ and sqrt() - all of which we have correctly implemented. I therefore believe that, for the remaining missing functions, the Project would be best served by committing the best code that is currently available under a suitable license and cleaning it up over time (as was done for the current libm). -- Peter Jeremy pgpPVXxJTjV0R.pgp Description: PGP signature
Re: Use of C99 extra long double math functions after r236148
On Friday, July 13, 2012 7:41:00 am Peter Jeremy wrote: AFAIK, none of the relevant standards (POSIX, IEEE754) have any precision requirements for functions other than +-*/ and sqrt() - all of which we have correctly implemented. I therefore believe that, for the remaining missing functions, the Project would be best served by committing the best code that is currently available under a suitable license and cleaning it up over time (as was done for the current libm). I concur. However, are there any patches that we can commit now? It seems that there are some things that could go into the tree now that will certainly reduce the concerns of the R folk. -- John Baldwin ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 13 Jul 2012, at 13:18, John Baldwin wrote: On Friday, July 13, 2012 7:41:00 am Peter Jeremy wrote: AFAIK, none of the relevant standards (POSIX, IEEE754) have any precision requirements for functions other than +-*/ and sqrt() - all of which we have correctly implemented. I therefore believe that, for the remaining missing functions, the Project would be best served by committing the best code that is currently available under a suitable license and cleaning it up over time (as was done for the current libm). I concur. As do I. I'd also point out that the ONLY requirement for long double according to the standard is that it has at least the same precision as double. Therefore, any implementation of these functions that is no worse that the double version is compliant. Once we have something meeting a minimum standard, then I'm very happy to see it improved, but having C99 functions missing now is just embarrassing while we're working on adding C11 features. David P.S. Someone said earlier that our clang still lacks some C99 features. Please point me at the relevant clang PRs and I'll be happy to work on them. There are quite a few open issues for C11 support, but C99 is, as far as I know, done. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Fri, Jul 13, 2012 at 09:41:00PM +1000, Peter Jeremy wrote: On 2012-Jul-11 15:32:47 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: I know an approach to implementing many of the missing functions. Are you willing to share this insight so someone else could do the work? For the missing trig and hyperbolic functions of complex float and complex double arguments, one can follow the approach in msun/src/s_ccosh[f].c. For the missing hyperbolic functions of long double and complex long double argument, people will need to wait for expl() and expm1l() to show up in tree. I have a working expl() and almost working expm1l() for ld80 systems. Until a few months ago, I did not have access to a ld128 system, so I haven't wasted my time implementing ld128 version that I could not test. When I do find some free time, I look at what is missing and start to put together a new function. At the moment, it seems that it takes 3+ years to get a new function written, tested, and committed. And, from what I can see, much of this is done quietly - which opens up the possibility that two people might both implement the same code or that people will avoid the area in fear of treading on someone else's toes. As I said previously, I believe the existing wiki page could be improved to form a central co-ordinating point to show what what activity is (or isn't) occurring. Well, I've posted the code I've written to freebsd-standards mailinglist more than once, and have only ever received comments from exactly 2 people. but most people seem to push the easy button and want to grab either cephes or netlib's libm. There are technical issues with this approach that I won't rehash again. Doing it properly requires significant effort by people with fairly specialised skills. Whilst the project has several people with the skills, it appears that none of them currently have the time. In the meantime, FreeBSD is taking free kicks from other FOSS groups that have gone down the quick-and-dirty path. AFAIK, none of the relevant standards (POSIX, IEEE754) have any precision requirements for functions other than +-*/ and sqrt() - all of which we have correctly implemented. I therefore believe that, for the remaining missing functions, the Project would be best served by committing the best code that is currently available under a suitable license and cleaning it up over time (as was done for the current libm). Well, hopefully, the someone who takes the best available code also tests this code before it is committed. I suspect that some if not all of those codes that involve complex arguments will have problems with the requirements in n1256.pdf[1], because neither gcc in base nor clang do complex arithmetic correctly under certains conditions when an expression uses I from complex.h. [1] http://www.open-std.org/jtc1/sc22/wg14/www/docs/n1256.pdf -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Fri, Jul 13, 2012 at 01:53:39PM +0100, David Chisnall wrote: On 13 Jul 2012, at 13:18, John Baldwin wrote: On Friday, July 13, 2012 7:41:00 am Peter Jeremy wrote: AFAIK, none of the relevant standards (POSIX, IEEE754) have any precision requirements for functions other than +-*/ and sqrt() - all of which we have correctly implemented. I therefore believe that, for the remaining missing functions, the Project would be best served by committing the best code that is currently available under a suitable license and cleaning it up over time (as was done for the current libm). I concur. As do I. I'd also point out that the ONLY requirement for long double according to the standard is that it has at least the same precision as double. Therefore, any implementation of these functions that is no worse that the double version is compliant. Once we have something meeting a minimum standard, then I'm very happy to see it improved, but having C99 functions missing now is just embarrassing while we're working on adding C11 features. I'd be curious how well the GPL functions in Linux compare to the NetBSD functions. I don't suppose we could grab some of the public domain routines in NetLib? David P.S. Someone said earlier that our clang still lacks some C99 features. Please point me at the relevant clang PRs and I'll be happy to work on them. There are quite a few open issues for C11 support, but C99 is, as far as I know, done. Diane -- - d...@freebsd.org d...@db.net http://www.db.net/~db Nowadays tar can compress using yesterdays latest technologies! ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Fri, Jul 13, 2012, David Chisnall wrote: As do I. I'd also point out that the ONLY requirement for long double according to the standard is that it has at least the same precision as double. Therefore, any implementation of these functions that is no worse that the double version is compliant. Once we have something meeting a minimum standard, then I'm very happy to see it improved, but having C99 functions missing now is just embarrassing while we're working on adding C11 features. There are several things wrong with this reasoning, but pragmatically the conclusion may be right: we do have a long list of users who would prefer a dubious implementation to none at all. I propose we set a timeframe for this, on the order of a few months. A rough outline might be something like: mid-August: expl logl log2l log10l -- just need to clean up Bruce and Steve's work; Steve recently sent me patches for expl, which I hope get committed soon mid-September: acoshl asinhl atanhl coshl sinhl tanhl -- easy once expl is in; others could probably help mid-October: powl expm1l mid-November: most complex.h functions If the schedule can't be met, then we can just import Cephes as an interim solution without further ado. This provides Bruce and Steve an opportunity to commit what they have been working on, without forcing the rest of the FreeBSD community to wait indefinitely for the pie in the sky. By the way, the trig and complex functions are areas where anyone with some calculus background could contribute. If anyone is interested in helping out, I'd be happy to coordinate things and review patches, although I will be unavailable for much of August. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 07/13/12 10:58, David Schultz wrote: On Fri, Jul 13, 2012, David Chisnall wrote: As do I. I'd also point out that the ONLY requirement for long double according to the standard is that it has at least the same precision as double. Therefore, any implementation of these functions that is no worse that the double version is compliant. Once we have something meeting a minimum standard, then I'm very happy to see it improved, but having C99 functions missing now is just embarrassing while we're working on adding C11 features. There are several things wrong with this reasoning, but pragmatically the conclusion may be right: we do have a long list of users who would prefer a dubious implementation to none at all. I propose we set a timeframe for this, on the order of a few months. A rough outline might be something like: mid-August: expl logl log2l log10l -- just need to clean up Bruce and Steve's work; Steve recently sent me patches for expl, which I hope get committed soon mid-September: acoshl asinhl atanhl coshl sinhl tanhl -- easy once expl is in; others could probably help mid-October: powl expm1l mid-November: most complex.h functions If the schedule can't be met, then we can just import Cephes as an interim solution without further ado. This provides Bruce and Steve an opportunity to commit what they have been working on, without forcing the rest of the FreeBSD community to wait indefinitely for the pie in the sky. This sounds fantastic. By the way, the trig and complex functions are areas where anyone with some calculus background could contribute. If anyone is interested in helping out, I'd be happy to coordinate things and review patches, although I will be unavailable for much of August. I would be happy to help. Just give me a sense of direction of what I should try and work on, when and as you are ready. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 13 July 2012 09:07, Stephen Montgomery-Smith step...@missouri.edu wrote: On 07/13/12 10:58, David Schultz wrote: On Fri, Jul 13, 2012, David Chisnall wrote: As do I. I'd also point out that the ONLY requirement for long double according to the standard is that it has at least the same precision as double. Therefore, any implementation of these functions that is no worse that the double version is compliant. Once we have something meeting a minimum standard, then I'm very happy to see it improved, but having C99 functions missing now is just embarrassing while we're working on adding C11 features. There are several things wrong with this reasoning, but pragmatically the conclusion may be right: we do have a long list of users who would prefer a dubious implementation to none at all. I propose we set a timeframe for this, on the order of a few months. A rough outline might be something like: mid-August: expl logl log2l log10l -- just need to clean up Bruce and Steve's work; Steve recently sent me patches for expl, which I hope get committed soon mid-September: acoshl asinhl atanhl coshl sinhl tanhl -- easy once expl is in; others could probably help mid-October: powl expm1l mid-November: most complex.h functions If the schedule can't be met, then we can just import Cephes as an interim solution without further ado. This provides Bruce and Steve an opportunity to commit what they have been working on, without forcing the rest of the FreeBSD community to wait indefinitely for the pie in the sky. +1 If we do import Cephes the questionable functions should probably be explicitly marked somewhere so that if there is still $someone can still work on them though. -- Eitan Adler ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
Just to jump back into the fray a bit, since this point hasn't been articulated well. On Jul 10, 2012, at 6:55 PM, Peter Jeremy wrote: On 2012-Jul-08 19:01:07 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: Well, on the most popular hardware (that being i386/amd64), ld80 will use hardware fp instruction while ld128 must be done completely in software. The speed difference is significant. AFAIK, of the architectures that FreeBSD supports, only sparc64 defines ld128 in the architecture and I don't believe there are any SPARC chip implementations that implement ld128 math in hardware. We shouldn't be gating the new math on an issue that only affects sparc64 machines. If they have ld80 level of support for that architecture, then that is sufficient to get things into the tree. There's no real benefit from making numerics good on sparc64 for the project, since our support for the platform isn't stellar and the platform itself is getting a bit long in the tooth. That said, if people want to do it, be my guest. If it is important enough to catch someone's attention, then it is important enough to have. It just isn't important enough to be a gating factor if nobody has signed up for it yet. Warner ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Friday, July 13, 2012 11:58:05 am David Schultz wrote: On Fri, Jul 13, 2012, David Chisnall wrote: As do I. I'd also point out that the ONLY requirement for long double according to the standard is that it has at least the same precision as double. Therefore, any implementation of these functions that is no worse that the double version is compliant. Once we have something meeting a minimum standard, then I'm very happy to see it improved, but having C99 functions missing now is just embarrassing while we're working on adding C11 features. There are several things wrong with this reasoning, but pragmatically the conclusion may be right: we do have a long list of users who would prefer a dubious implementation to none at all. I propose we set a timeframe for this, on the order of a few months. A rough outline might be something like: mid-August: expl logl log2l log10l -- just need to clean up Bruce and Steve's work; Steve recently sent me patches for expl, which I hope get committed soon mid-September: acoshl asinhl atanhl coshl sinhl tanhl -- easy once expl is in; others could probably help mid-October: powl expm1l mid-November: most complex.h functions If the schedule can't be met, then we can just import Cephes as an interim solution without further ado. This provides Bruce and Steve an opportunity to commit what they have been working on, without forcing the rest of the FreeBSD community to wait indefinitely for the pie in the sky. By the way, the trig and complex functions are areas where anyone with some calculus background could contribute. If anyone is interested in helping out, I'd be happy to coordinate things and review patches, although I will be unavailable for much of August. I think this sounds like an excellent plan, thanks! -- John Baldwin ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 2012-Jul-13 11:58:05 -0400, David Schultz d...@freebsd.org wrote: I propose we set a timeframe for this, on the order of a few months. ... If the schedule can't be met, then we can just import Cephes as an interim solution without further ado. This provides Bruce and Steve an opportunity to commit what they have been working on, without forcing the rest of the FreeBSD community to wait indefinitely for the pie in the sky. This sounds good to me as well and I'd be happy to help. -- Peter Jeremy pgpmY7CNvs676.pgp Description: PGP signature
Re: Use of C99 extra long double math functions after r236148
On Jul 13, 2012, at 4:38 PM, Peter Jeremy wrote: On 2012-Jul-13 11:58:05 -0400, David Schultz d...@freebsd.org wrote: I propose we set a timeframe for this, on the order of a few months. ... If the schedule can't be met, then we can just import Cephes as an interim solution without further ado. This provides Bruce and Steve an opportunity to commit what they have been working on, without forcing the rest of the FreeBSD community to wait indefinitely for the pie in the sky. This sounds good to me as well and I'd be happy to help. Also sounds good to me. If Peter gets busy, I'd be able to help as well with the back-stop proposal. Warner ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Sat, Jul 14, 2012 at 08:38:08AM +1000, Peter Jeremy wrote: On 2012-Jul-13 11:58:05 -0400, David Schultz d...@freebsd.org wrote: I propose we set a timeframe for this, on the order of a few months. ... If the schedule can't be met, then we can just import Cephes as an interim solution without further ado. This provides Bruce and Steve an opportunity to commit what they have been working on, without forcing the rest of the FreeBSD community to wait indefinitely for the pie in the sky. This sounds good to me as well and I'd be happy to help. Perfect, all that is needed. I'd also be happy to help of course. -- Peter Jeremy Diane -- - d...@freebsd.org d...@db.net http://www.db.net/~db Nowadays tar can compress using yesterdays latest technologies! ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 14/07/2012, at 2:35, Eitan Adler wrote: If the schedule can't be met, then we can just import Cephes as an interim solution without further ado. This provides Bruce and Steve an opportunity to commit what they have been working on, without forcing the rest of the FreeBSD community to wait indefinitely for the pie in the sky. +1 If we do import Cephes the questionable functions should probably be explicitly marked somewhere so that if there is still $someone can still work on them though. Perhaps the same way mktemp et al are marked.. -- Daniel O'Connor software and network engineer for Genesis Software - http://www.gsoft.com.au The nice thing about standards is that there are so many of them to choose from. -- Andrew Tanenbaum GPG Fingerprint - 5596 B766 97C0 0E94 4347 295E E593 DC20 7B3F CE8C
Re: Use of C99 extra long double math functions after r236148
On Tue, Jul 10, 2012 at 08:02:33AM -0600, Warner Losh wrote: On Jul 10, 2012, at 3:10 AM, Rainer Hurling wrote: Not having R would be a bit pain in my backside. That's one of the practical considerations that I was talking about. It is very real, and if I have to, I'll commit the #define junk I railed against to get it back. Please, let's get some progress. I have some time to help. Not having complex functions is a PITA for software defined radio programs... I submitted this PR in frustration. I know, I meant to close it already. Number: 147599 Category: kern Synopsis: [libm] [patch] Import netbsd complex functions into our libm Warner Diane -- - d...@freebsd.org d...@db.net http://www.db.net/~db Nowadays tar can compress using yesterdays latest technologies! ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Wed, Jul 11, 2012 at 04:20:09PM -0500, Diane Bruce wrote: On Tue, Jul 10, 2012 at 08:02:33AM -0600, Warner Losh wrote: On Jul 10, 2012, at 3:10 AM, Rainer Hurling wrote: Not having R would be a bit pain in my backside. That's one of the practical considerations that I was talking about. It is very real, and if I have to, I'll commit the #define junk I railed against to get it back. Please, let's get some progress. I have some time to help. Not having complex functions is a PITA for software defined radio programs... I submitted this PR in frustration. I know, I meant to close it already. Number: 147599 Category: kern Synopsis: [libm] [patch] Import netbsd complex functions into our libm You submitted on June 6th, 2010. I commented on why the patch should be avoided on June 15th, 2010. I see no follow-ups from you that give the details on how you went about testing and fixing the code? -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Wed, Jul 11, 2012 at 02:43:46PM -0700, Steve Kargl wrote: On Wed, Jul 11, 2012 at 04:20:09PM -0500, Diane Bruce wrote: On Tue, Jul 10, 2012 at 08:02:33AM -0600, Warner Losh wrote: On Jul 10, 2012, at 3:10 AM, Rainer Hurling wrote: ... You submitted on June 6th, 2010. I commented on why the patch should be avoided on June 15th, 2010. I see no follow-ups from you that give the details on how you went about testing and fixing the code? Steve, I misunderstood. Someone told me you and Bruce had some code which made this PR unnecessary, I did mean to close the PR as this is what I had heard. If that is not true I'd be happy to help you guys in any way I can. And my apologies, I should have verfied directly with you rather than listening to hearsay. -- Steve Diane -- - d...@freebsd.org d...@db.net http://www.db.net/~db Nowadays tar can compress using yesterdays latest technologies! ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Wed, Jul 11, 2012 at 04:54:14PM -0500, Diane Bruce wrote: On Wed, Jul 11, 2012 at 02:43:46PM -0700, Steve Kargl wrote: On Wed, Jul 11, 2012 at 04:20:09PM -0500, Diane Bruce wrote: On Tue, Jul 10, 2012 at 08:02:33AM -0600, Warner Losh wrote: On Jul 10, 2012, at 3:10 AM, Rainer Hurling wrote: ... You submitted on June 6th, 2010. I commented on why the patch should be avoided on June 15th, 2010. I see no follow-ups from you that give the details on how you went about testing and fixing the code? Steve, I misunderstood. Someone told me you and Bruce had some code which made this PR unnecessary, I did mean to close the PR as this is what I had heard. If that is not true I'd be happy to help you guys in any way I can. And my apologies, I should have verfied directly with you rather than listening to hearsay. I know an approach to implementing many of the missing functions. The problem is that over the last 2 years or so I have had very little time to work on coding nor am I paid to work on such code. When I do find some free time, I look at what is missing and start to put together a new function. At the moment, it seems that it takes 3+ years to get a new function written, tested, and committed. Given that it seems that only das@, bde@, and I work on libm, I suspect it will be completed by the time I retire. At one time I had hoped others would step up to help write the missing function, but most people seem to push the easy button and want to grab either cephes or netlib's libm. There are technical issues with this approach that I won't rehash again. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 09.07.2012 01:36 (UTC+2), Steve Kargl wrote: On Sun, Jul 08, 2012 at 11:51:56AM -0600, Warner Losh wrote: On Jul 8, 2012, at 6:40 AM, David Schultz wrote: On Tue, May 29, 2012, Peter Jeremy wrote: On 2012-May-28 15:54:06 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: Given that cephes was written years before C99 was even conceived, I suspect all functions are sub-standard. Well, most of cephes was written before C99. The C99 parts of cephes were written to turn it into a complete C99 implementation. I'm a bit late to the party, but I thought I'd chime in with some context. We did consider using Cephes years ago, and even got permission from the author to release it under an acceptable license. We later decided not to use it for technical reasons. By the way, virtually none of the people who have complained about the missing functions actually need them. Mostly they just want to compile some software that was written by a naive programmer who thought it would be cool to use the most precise type available. The complex functions are even less commonly needed, and the truth is that they have no business being part of the C standard anyway. The question remains of what to do about the missing functions. Bruce and Steve have been working on expl and logl for years. If those ever get in the tree, the remaining long double functions are easy. Those functions are basically done, modulo a bunch of cleanup and testing, and I encourage any mathematically inclined folks who are interested in pushing things along to get in touch with them. I'm not going to have any time myself for a few months at least. Where can I find these? I've posted expl() a few times for the ld80 version. I don't have an ld128 version, which is why I have yet to submit a formal patch for expl(). I also have an ld80 expm1l(). I have a copy of bde's ld80 logl(). IIRC, bde wrote an ld128, but I don't have nor do I know if it has been tested. Do you think your version from http://www.freebsd.org/cgi/query-pr.cgi?pr=152415 for expl() ld80 version could be the one getting into head? Would you be willing to commit it? As far as I understand from discussions on R mailing list (r-de...@r-project.org), they plan to reduce the emulation and/or workaround of long and complex math functions for FreeBSD and other systems with their next releases of R devel. So we could really need some progress with our C99 conform math functions ;-) Lastly, there's the question of mediocre alternatives, such as solutions that get the boundary cases wrong or don't handle 128-bit floating point. For the exponential and logarithmic functions, Bruce and Steve have already written good implementations, so there's no reason to settle for less. As for the other long double functions, bringing in some Cephes code in a separate directory as a temporary fix might be the way to go. I don't like that solution, and Steve raises some good technical points about why it isn't ideal; however, a better solution is more than a decade overdue, and people are justified in finding that unacceptable. Don't let the perfect be the enemy of the good. It is better to have OK implementations of these functions than none at all. You'll need to define OK, because some of the implementations I've read are poor. I also hope that before anything is committed to libm that the person doing the committing does some rather thorough testing of the code. We originally had so-so double support, but bruce spent many years optimizing them to make them very good. If we were just starting out, and hadn't let 10 years get behind us, I'd give the substandard argument some weight. But now that we're 13 years down the line from c99's publication I think we need to relax our standards a bit. I'd even argue that these functions being a little bad could easily spur people to make them better. Their absence makes people just #define llexp(x) lexp(x), etc. :( Of course, the converse argument is that people have had 10 years to learn enough about floating point to actually write and submit code. I knew very little about FP before I wrote sqrtl(). I had a need for sqrtl() and so I went looking for a solution. PS: I also wrote sincos[fl](), which is very handy for the complex trig functions. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Tue, Jul 10, 2012 at 11:10:05AM +0200, Rainer Hurling wrote: Do you think your version from http://www.freebsd.org/cgi/query-pr.cgi?pr=152415 for expl() ld80 version could be the one getting into head? Would you be willing to commit it? That's a fairly early version of the ld80 expl(). bde and das reviewed the code and made several suggested changes. As far as I understand from discussions on R mailing list (r-de...@r-project.org), they plan to reduce the emulation and/or workaround of long and complex math functions for FreeBSD and other systems with their next releases of R devel. So we could really need some progress with our C99 conform math functions ;-) I don't know R, and I don't understand what is meant by 'they plan to reduce the emulation and/or workaround of ... funtions'. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Jul 10, 2012, at 3:10 AM, Rainer Hurling wrote: As far as I understand from discussions on R mailing list (r-de...@r-project.org), they plan to reduce the emulation and/or workaround of long and complex math functions for FreeBSD and other systems with their next releases of R devel. So we could really need some progress with our C99 conform math functions ;-) Not having R would be a bit pain in my backside. That's one of the practical considerations that I was talking about. It is very real, and if I have to, I'll commit the #define junk I railed against to get it back. Please, let's get some progress. I have some time to help. Warner ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Jul 10, 2012, at 7:45 AM, Steve Kargl wrote: On Tue, Jul 10, 2012 at 11:10:05AM +0200, Rainer Hurling wrote: Do you think your version from http://www.freebsd.org/cgi/query-pr.cgi?pr=152415 for expl() ld80 version could be the one getting into head? Would you be willing to commit it? That's a fairly early version of the ld80 expl(). bde and das reviewed the code and made several suggested changes. Cool. Send them to me and I'll move them into the tree. As far as I understand from discussions on R mailing list (r-de...@r-project.org), they plan to reduce the emulation and/or workaround of long and complex math functions for FreeBSD and other systems with their next releases of R devel. So we could really need some progress with our C99 conform math functions ;-) I don't know R, and I don't understand what is meant by 'they plan to reduce the emulation and/or workaround of ... funtions'. R is a totally awesome statistical and graphics packaging. It means they are dropping support for freeBSD because we don't have these functions. That totally sucks, and to prevent that, I'd be willing to commit the #define nonsense. Unless somebody has something better they could forward me for integration... Warner___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 10.07.2012 16:02 (UTC+2), Warner Losh wrote: On Jul 10, 2012, at 3:10 AM, Rainer Hurling wrote: As far as I understand from discussions on R mailing list (r-de...@r-project.org), they plan to reduce the emulation and/or workaround of long and complex math functions for FreeBSD and other systems with their next releases of R devel. So we could really need some progress with our C99 conform math functions ;-) Not having R would be a bit pain in my backside. That's one of the practical considerations that I was talking about. It is very real, and if I have to, I'll commit the #define junk I railed against to get it back. Please, let's get some progress. I have some time to help. Yes, thank you Warner, that is also my problem. As I wrote some weeks ago (05/28/2012) when starting this thread, I am using FreeBSD as a scientific desktop because of its good scaling properties. For some years now, FreeBSD fits all our needs with R, SAGA GIS, PostgreSQL and some more. If I would not be able to run upcoming versions of R on FreeBSD any more, that would be really, really hard :-( Rainer Warner ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Tue, Jul 10, 2012, Rainer Hurling wrote: On 10.07.2012 16:02 (UTC+2), Warner Losh wrote: On Jul 10, 2012, at 3:10 AM, Rainer Hurling wrote: As far as I understand from discussions on R mailing list (r-de...@r-project.org), they plan to reduce the emulation and/or workaround of long and complex math functions for FreeBSD and other systems with their next releases of R devel. So we could really need some progress with our C99 conform math functions ;-) Not having R would be a bit pain in my backside. That's one of the practical considerations that I was talking about. It is very real, and if I have to, I'll commit the #define junk I railed against to get it back. Please, let's get some progress. I have some time to help. Yes, thank you Warner, that is also my problem. As I wrote some weeks ago (05/28/2012) when starting this thread, I am using FreeBSD as a scientific desktop because of its good scaling properties. For some years now, FreeBSD fits all our needs with R, SAGA GIS, PostgreSQL and some more. If I would not be able to run upcoming versions of R on FreeBSD any more, that would be really, really hard :-( Do you have a list of the essential functions here? There are 17 long double functions and some complex functions missing, but only a handful of those are of general interest. The reason I ask is that if R is just looking for a few missing functions that are already mostly implemented, then the best solution is probably to finish that work. But if it's expecting us to have something arcane like long double Bessel functions of the first kind, then we need to pursue a workaround in the short term. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Jul 10, 2012, at 9:11 AM, David Schultz wrote: On Tue, Jul 10, 2012, Rainer Hurling wrote: On 10.07.2012 16:02 (UTC+2), Warner Losh wrote: On Jul 10, 2012, at 3:10 AM, Rainer Hurling wrote: As far as I understand from discussions on R mailing list (r-de...@r-project.org), they plan to reduce the emulation and/or workaround of long and complex math functions for FreeBSD and other systems with their next releases of R devel. So we could really need some progress with our C99 conform math functions ;-) Not having R would be a bit pain in my backside. That's one of the practical considerations that I was talking about. It is very real, and if I have to, I'll commit the #define junk I railed against to get it back. Please, let's get some progress. I have some time to help. Yes, thank you Warner, that is also my problem. As I wrote some weeks ago (05/28/2012) when starting this thread, I am using FreeBSD as a scientific desktop because of its good scaling properties. For some years now, FreeBSD fits all our needs with R, SAGA GIS, PostgreSQL and some more. If I would not be able to run upcoming versions of R on FreeBSD any more, that would be really, really hard :-( Do you have a list of the essential functions here? There are 17 long double functions and some complex functions missing, but only a handful of those are of general interest. The reason I ask is that if R is just looking for a few missing functions that are already mostly implemented, then the best solution is probably to finish that work. But if it's expecting us to have something arcane like long double Bessel functions of the first kind, then we need to pursue a workaround in the short term. I'll get a list, but they are things like long long double exp and ln. I don't think they were the Bessel functions. Warner ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 07/09/2012 12:02 AM, Steve Kargl wrote: Yep. Another example is the use of upward recurion to compute Bessel functions where the argument is larger than the order. The algorithm is known to be unstable. By upward recursion, do you mean equation (1) in http://mathworld.wolfram.com/BesselFunction.html? So what do people use. Maybe something like http://en.wikipedia.org/wiki/Bessel_function#Asymptotic_forms (second set of equations), but finding some asymptotics with a few extra terms in them? ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Tue, Jul 10, 2012 at 11:39:59AM -0500, Stephen Montgomery-Smith wrote: On 07/09/2012 12:02 AM, Steve Kargl wrote: Yep. Another example is the use of upward recurion to compute Bessel functions where the argument is larger than the order. The algorithm is known to be unstable. By upward recursion, do you mean equation (1) in http://mathworld.wolfram.com/BesselFunction.html? Yes. So what do people use. Maybe something like http://en.wikipedia.org/wiki/Bessel_function#Asymptotic_forms (second set of equations), but finding some asymptotics with a few extra terms in them? They use downward recursion, which is known to be stable. NIST has revised Abramowitz and Stegun, and it is available on line. For Bessel function computations, look at http://dlmf.nist.gov/10.74 and more importantly example 1 under the following link http://dlmf.nist.gov/3.6#v -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 10.07.2012 17:11 (UTC+2), David Schultz wrote: On Tue, Jul 10, 2012, Rainer Hurling wrote: On 10.07.2012 16:02 (UTC+2), Warner Losh wrote: On Jul 10, 2012, at 3:10 AM, Rainer Hurling wrote: As far as I understand from discussions on R mailing list (r-de...@r-project.org), they plan to reduce the emulation and/or workaround of long and complex math functions for FreeBSD and other systems with their next releases of R devel. So we could really need some progress with our C99 conform math functions ;-) Not having R would be a bit pain in my backside. That's one of the practical considerations that I was talking about. It is very real, and if I have to, I'll commit the #define junk I railed against to get it back. Please, let's get some progress. I have some time to help. Yes, thank you Warner, that is also my problem. As I wrote some weeks ago (05/28/2012) when starting this thread, I am using FreeBSD as a scientific desktop because of its good scaling properties. For some years now, FreeBSD fits all our needs with R, SAGA GIS, PostgreSQL and some more. If I would not be able to run upcoming versions of R on FreeBSD any more, that would be really, really hard :-( Do you have a list of the essential functions here? There are 17 long double functions and some complex functions missing, but only a handful of those are of general interest. The reason I ask is that if R is just looking for a few missing functions that are already mostly implemented, then the best solution is probably to finish that work. But if it's expecting us to have something arcane like long double Bessel functions of the first kind, then we need to pursue a workaround in the short term. That is, what I found by grepping the sources of a recent R development version: expl: src/nmath/pnchisq.c log10l: src/extra/trio/trio.c log1pl: src/nmath/pnbeta.c logl: src/nmath/dnbeta.c src/nmath/pnbeta.c powl: src/extra/trio/triostr.c src/extra/trio/trio.c src/main/format.c In the R NEWS file I found the following three entries about C99 (complex) functions: NEWS:l2044 The C99 functions acosh, asinh, atanh, snprintf and vsnprintf are now required. NEWS:l3032 The C99 double complex type is now required. The C99 complex trigonometric functions (such as csin) are not currently required (FreeBSD lacks most of them): substitutes are used if they are missing. NEWS:l3277 Complex arithmetic (notably z^n for complex z and integer n) gave incorrect results since R 2.10.0 on platforms without C99 complex support. This and some lesser issues in trigonometric functions have been corrected. Such platforms were rare (we know of Cygwin and FreeBSD). However, because of new compiler optimizations in the way complex arguments are handled, the same code was selected on x86_64 Linux with gcc 4.5.x at the default -O2 optimization (but not at -O). From the discussions on the R-devel mailing list we could expect, that there will be more complex (long, longlong) functions introduced in the near future. (But it is very likely that b.f. does know much more about this.) Rainer BTW: There seems to be a discrepancy about missing functions listed in http://wiki.freebsd.org/MissingMathStuff and in http://svnweb.freebsd.org/base/head/lib/msun/src/math.h?r1=227472r2=236148pathrev=236148. So the wiki is a bit outdated now? ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 07/10/2012 11:50 AM, Steve Kargl wrote: On Tue, Jul 10, 2012 at 11:39:59AM -0500, Stephen Montgomery-Smith wrote: On 07/09/2012 12:02 AM, Steve Kargl wrote: Yep. Another example is the use of upward recurion to compute Bessel functions where the argument is larger than the order. The algorithm is known to be unstable. By upward recursion, do you mean equation (1) in http://mathworld.wolfram.com/BesselFunction.html? Yes. So what do people use. Maybe something like http://en.wikipedia.org/wiki/Bessel_function#Asymptotic_forms (second set of equations), but finding some asymptotics with a few extra terms in them? They use downward recursion, which is known to be stable. NIST has revised Abramowitz and Stegun, and it is available on line. For Bessel function computations, look at http://dlmf.nist.gov/10.74 and more importantly example 1 under the following link http://dlmf.nist.gov/3.6#v These algorithms are going to be subject to the same problems as using Taylor's series to evaluate exp(x) for x0. The computation will require several floating point operations. Even if the method is not unstable, I would think getting a ULP of about 1 is next to impossible, that is, unless one performs all the computations in a higher precision, and then rounds at the end. Whereas as a scientist, having a bessel function computed to within 10 ULP would be just fine. I am looking at the man page of j0 for FreeBSD-8.3. It says in conforms to IEEE Std 1003.1-2001. Does that specify a certain ULP? I am searching around in this document, and I am finding nothing. Whereas the IEEE-754 document seems rather rigid, but on the other hand it doesn't specifically talk about math functions other than sqrt. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Tue, Jul 10, 2012 at 01:11:38PM -0500, Stephen Montgomery-Smith wrote: On 07/10/2012 11:50 AM, Steve Kargl wrote: On Tue, Jul 10, 2012 at 11:39:59AM -0500, Stephen Montgomery-Smith wrote: On 07/09/2012 12:02 AM, Steve Kargl wrote: Yep. Another example is the use of upward recurion to compute Bessel functions where the argument is larger than the order. The algorithm is known to be unstable. By upward recursion, do you mean equation (1) in http://mathworld.wolfram.com/BesselFunction.html? Yes. So what do people use. Maybe something like http://en.wikipedia.org/wiki/Bessel_function#Asymptotic_forms (second set of equations), but finding some asymptotics with a few extra terms in them? They use downward recursion, which is known to be stable. NIST has revised Abramowitz and Stegun, and it is available on line. For Bessel function computations, look at http://dlmf.nist.gov/10.74 and more importantly example 1 under the following link http://dlmf.nist.gov/3.6#v These algorithms are going to be subject to the same problems as using Taylor's series to evaluate exp(x) for x0. The computation will require several floating point operations. Even if the method is not unstable, I would think getting a ULP of about 1 is next to impossible, that is, unless one performs all the computations in a higher precision, and then rounds at the end. Downward recurion is used for n 1. j0(x) and j1(x) are computed using two different schemes. If x = 2, then the ascending series can be summed (this only takes somewhere around 10 terms for double). Alternately, one can use a minimax polynomial for x in [0:2]. For x 2, the Hankel expansions are used (see 10.17.3 at http://dlmf.nist.gov/10.17). The summations are approximated by minimax polynomials. See msun/src/s_j0.c. I've tested j0f() and know that one can get less than 1 ULP over various ranges of x. I also know that when is close to a zero of j0(x), then one has ULP on the order of 100. Whereas as a scientist, having a bessel function computed to within 10 ULP would be just fine. As someone who has spent a portion of career computing Bessel functions, I am disinclined to agree with your statement. I am looking at the man page of j0 for FreeBSD-8.3. It says in conforms to IEEE Std 1003.1-2001. Does that specify a certain ULP? I doubt that it would specify ULP for any function other than perhaps sqrt(). I don't have IEEE Std 1003.1-2001 handy, but I would guess that it simply statements the math funtion return an floating point approximation to its true value. I am searching around in this document, and I am finding nothing. Whereas the IEEE-754 document seems rather rigid, but on the other hand it doesn't specifically talk about math functions other than sqrt. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Tue, Jul 10, 2012, Rainer Hurling wrote: On 10.07.2012 17:11 (UTC+2), David Schultz wrote: On Tue, Jul 10, 2012, Rainer Hurling wrote: On 10.07.2012 16:02 (UTC+2), Warner Losh wrote: On Jul 10, 2012, at 3:10 AM, Rainer Hurling wrote: As far as I understand from discussions on R mailing list (r-de...@r-project.org), they plan to reduce the emulation and/or workaround of long and complex math functions for FreeBSD and other systems with their next releases of R devel. So we could really need some progress with our C99 conform math functions ;-) Not having R would be a bit pain in my backside. That's one of the practical considerations that I was talking about. It is very real, and if I have to, I'll commit the #define junk I railed against to get it back. Please, let's get some progress. I have some time to help. Yes, thank you Warner, that is also my problem. As I wrote some weeks ago (05/28/2012) when starting this thread, I am using FreeBSD as a scientific desktop because of its good scaling properties. For some years now, FreeBSD fits all our needs with R, SAGA GIS, PostgreSQL and some more. If I would not be able to run upcoming versions of R on FreeBSD any more, that would be really, really hard :-( Do you have a list of the essential functions here? There are 17 long double functions and some complex functions missing, but only a handful of those are of general interest. The reason I ask is that if R is just looking for a few missing functions that are already mostly implemented, then the best solution is probably to finish that work. But if it's expecting us to have something arcane like long double Bessel functions of the first kind, then we need to pursue a workaround in the short term. That is, what I found by grepping the sources of a recent R development version: expl: src/nmath/pnchisq.c logl: src/nmath/dnbeta.c src/nmath/pnbeta.c Bruce has versions of these that could be committed with some cleanup. It's a matter of sorting through about 1200 emails from him and 3 source trees to find the most up-to-date patches, then cleaning them up and testing and committing them. I have no time right now, but I will do at least the first step as soon as I can, and try to get the patches to someone willing to do the final few steps. log10l: src/extra/trio/trio.c log1pl: src/nmath/pnbeta.c If Bruce doesn't already have implementations of these, they are easy wrappers around logl() or some internal k_logl in Bruce's implementation. powl: src/extra/trio/triostr.c src/extra/trio/trio.c src/main/format.c It's hard to do a good job on powl(), but the simple approach (exp(log(x)*y)) plus a few special cases may suffice for many uses. NEWS:l2044 The C99 functions acosh, asinh, atanh, snprintf and vsnprintf are now required. We have had them forever. NEWS:l3032 The C99 double complex type is now required. The C99 complex trigonometric functions (such as csin) are not currently required (FreeBSD lacks most of them): substitutes are used if they are missing. We have these (but not the inverse trig functions). NEWS:l3277 Complex arithmetic (notably z^n for complex z and integer n) gave incorrect results since R 2.10.0 on platforms without C99 complex support. This and some lesser issues in trigonometric functions have been corrected. Such platforms were rare (we know of Cygwin and FreeBSD). However, because of new compiler optimizations in the way complex arguments are handled, the same code was selected on x86_64 Linux with gcc 4.5.x at the default -O2 optimization (but not at -O). Not sure if this is relevant. BTW: There seems to be a discrepancy about missing functions listed in http://wiki.freebsd.org/MissingMathStuff and in http://svnweb.freebsd.org/base/head/lib/msun/src/math.h?r1=227472r2=236148pathrev=236148. So the wiki is a bit outdated now? My list: REAL FUNCTIONS (17): long double log2l(long double); long double logl(long double); long double log1pl(long double); long double acoshl(long double); long double asinhl(long double); long double atanhl(long double); long double log10l(long double); long double expl(long double); long double expm1l(long double); long double coshl(long double); long double sinhl(long double); long double tanhl(long double); long double erfcl(long double); long double erfl(long double); long double powl(long double, long double); long double lgammal(long double); long double tgammal(long double); COMPLEX FUNCTIONS (37): long double complex cexpl(long double complex); long double complex ccosl(long double complex); long double complex ccoshl(long double complex); long double complex csinl(long double complex); long double complex csinhl(long double complex); long double complex ctanl(long double complex); long double complex ctanhl(long double
Re: Use of C99 extra long double math functions after r236148
On 2012-Jul-08 19:01:07 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: Well, on the most popular hardware (that being i386/amd64), ld80 will use hardware fp instruction while ld128 must be done completely in software. The speed difference is significant. AFAIK, of the architectures that FreeBSD supports, only sparc64 defines ld128 in the architecture and I don't believe there are any SPARC chip implementations that implement ld128 math in hardware. For that matter, I don't believe anything except x86 provides full IEEE FP support in hardware - most architectures require software assistance for subnormals and some corner cases. If your application happens to hit those cases often, performance will also suffer. On 2012-Jul-08 20:05:04 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: AFAIK, neither gcc in base nor clang would be c99 complaint even if all of the c99 math functions were available. That sort of argument can easily get circular. Lets get the C99 bits of libm out of the way and then we can have another bikeshed about the shortcomings of the compiler(s). On 2012-Jul-08 19:56:52 -0400, David Schultz d...@freebsd.org wrote: Yes, Bruce has ld128 versions, and clusteradm very kindly got us a sparc64 machine to test on. That was about the time I ran out of time to keep working on it. If someone wants to pick it up, that would be great. I have access to a couple of SPARC systems as well and would be willing to help work on the missing bits. On 2012-Jul-10 18:58:01 -0400, David Schultz d...@freebsd.org wrote: On Tue, Jul 10, 2012, Rainer Hurling wrote: powl: src/extra/trio/triostr.c src/extra/trio/trio.c src/main/format.c It's hard to do a good job on powl(), but the simple approach (exp(log(x)*y)) plus a few special cases may suffice for many uses. A simplistic exp(log(x)*y) throws away 15 bits of precision (size of the FP exponent field). cephes has a powl() that appears to do better or, alternatively, it shouldn't be too difficult to extend the approach used by __ieee754_pow() using long doubles. BTW: There seems to be a discrepancy about missing functions listed in http://wiki.freebsd.org/MissingMathStuff and in http://svnweb.freebsd.org/base/head/lib/msun/src/math.h?r1=227472r2=236148pathrev=236148. So the wiki is a bit outdated now? My list: [elided] I was thinking that a wiki page would be a good spot to co-ordinate the work (as well as making it clear what is still to be done). The existing page needs some TLC to be useful. -- Peter Jeremy pgpJMDQgZRF8K.pgp Description: PGP signature
Re: Use of C99 extra long double math functions after r236148
Hi, On 9 Jul 2012, at 06:02, Steve Kargl wrote: If you're doing accouting, hopefully, you're using BCD. Would be nice, but it's far too slow for financial analysis; the chip designers go out of their way to provide fast floating point. Fortunately 53 bits is usually plenty with 2dp max. -- Bob Bishop r...@gid.co.uk ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Sunday, July 08, 2012 7:56:52 pm David Schultz wrote: On Sun, Jul 08, 2012, Steve Kargl wrote: The question remains of what to do about the missing functions. Bruce and Steve have been working on expl and logl for years. If those ever get in the tree, the remaining long double functions are easy. Those functions are basically done, modulo a bunch of cleanup and testing, and I encourage any mathematically inclined folks who are interested in pushing things along to get in touch with them. I'm not going to have any time myself for a few months at least. Where can I find these? I've posted expl() a few times for the ld80 version. I don't have an ld128 version, which is why I have yet to submit a formal patch for expl(). I also have an ld80 expm1l(). I have a copy of bde's ld80 logl(). IIRC, bde wrote an ld128, but I don't have nor do I know if it has been tested. Yes, Bruce has ld128 versions, and clusteradm very kindly got us a sparc64 machine to test on. That was about the time I ran out of time to keep working on it. If someone wants to pick it up, that would be great. PS: I also wrote sincos[fl](), which is very handy for the complex trig functions. Yes, you should commit that! So it sounds like there are at least some outstanding patches that should go into the tree. Can we at least get these things committed (even if they are only enabled for certain platforms for now)? -- John Baldwin ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Tue, May 29, 2012, Peter Jeremy wrote: On 2012-May-28 15:54:06 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: Given that cephes was written years before C99 was even conceived, I suspect all functions are sub-standard. Well, most of cephes was written before C99. The C99 parts of cephes were written to turn it into a complete C99 implementation. I'm a bit late to the party, but I thought I'd chime in with some context. We did consider using Cephes years ago, and even got permission from the author to release it under an acceptable license. We later decided not to use it for technical reasons. By the way, virtually none of the people who have complained about the missing functions actually need them. Mostly they just want to compile some software that was written by a naive programmer who thought it would be cool to use the most precise type available. The complex functions are even less commonly needed, and the truth is that they have no business being part of the C standard anyway. The question remains of what to do about the missing functions. Bruce and Steve have been working on expl and logl for years. If those ever get in the tree, the remaining long double functions are easy. Those functions are basically done, modulo a bunch of cleanup and testing, and I encourage any mathematically inclined folks who are interested in pushing things along to get in touch with them. I'm not going to have any time myself for a few months at least. Lastly, there's the question of mediocre alternatives, such as solutions that get the boundary cases wrong or don't handle 128-bit floating point. For the exponential and logarithmic functions, Bruce and Steve have already written good implementations, so there's no reason to settle for less. As for the other long double functions, bringing in some Cephes code in a separate directory as a temporary fix might be the way to go. I don't like that solution, and Steve raises some good technical points about why it isn't ideal; however, a better solution is more than a decade overdue, and people are justified in finding that unacceptable. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Jul 8, 2012, at 6:40 AM, David Schultz wrote: On Tue, May 29, 2012, Peter Jeremy wrote: On 2012-May-28 15:54:06 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: Given that cephes was written years before C99 was even conceived, I suspect all functions are sub-standard. Well, most of cephes was written before C99. The C99 parts of cephes were written to turn it into a complete C99 implementation. I'm a bit late to the party, but I thought I'd chime in with some context. We did consider using Cephes years ago, and even got permission from the author to release it under an acceptable license. We later decided not to use it for technical reasons. By the way, virtually none of the people who have complained about the missing functions actually need them. Mostly they just want to compile some software that was written by a naive programmer who thought it would be cool to use the most precise type available. The complex functions are even less commonly needed, and the truth is that they have no business being part of the C standard anyway. The question remains of what to do about the missing functions. Bruce and Steve have been working on expl and logl for years. If those ever get in the tree, the remaining long double functions are easy. Those functions are basically done, modulo a bunch of cleanup and testing, and I encourage any mathematically inclined folks who are interested in pushing things along to get in touch with them. I'm not going to have any time myself for a few months at least. Where can I find these? Lastly, there's the question of mediocre alternatives, such as solutions that get the boundary cases wrong or don't handle 128-bit floating point. For the exponential and logarithmic functions, Bruce and Steve have already written good implementations, so there's no reason to settle for less. As for the other long double functions, bringing in some Cephes code in a separate directory as a temporary fix might be the way to go. I don't like that solution, and Steve raises some good technical points about why it isn't ideal; however, a better solution is more than a decade overdue, and people are justified in finding that unacceptable. Don't let the perfect be the enemy of the good. It is better to have OK implementations of these functions than none at all. We originally had so-so double support, but bruce spent many years optimizing them to make them very good. If we were just starting out, and hadn't let 10 years get behind us, I'd give the substandard argument some weight. But now that we're 13 years down the line from c99's publication I think we need to relax our standards a bit. I'd even argue that these functions being a little bad could easily spur people to make them better. Their absence makes people just #define llexp(x) lexp(x), etc. :( Warner Warner___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
Here is a technical question. I know that people always talk about ulp's in the context of how good a function implementation is. I think the ulp is the number of base 2 digits at the end of the mantissa that we cannot rely on. So if one were to write a naive implementation of lexp(x) that used Taylor's series if x is positive, and 1/lexp(-x) is x is negative - well one could fairly easily estimate an upper bound on the ulp, but it wouldn't be low (like ulp=1 or 2), but probably rather higher (ulp of the order of 10 or 20). So do people really work hard to get that last drop of ulp out of their calculations? Would a ulp=10 be considered unacceptable? Also, looking through the source code for the FreeBSD implementation of exp, I saw that they used some rather smart rational function. (I don't know how they came up with it.) Presumably a big part of the issue is to make the functions work rather fast. And a naive implementation of Taylor's series wouldn't be fast. But if people want lexp rather than exp, they must have already decided that accuracy is more important than speed. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Sun, Jul 08, 2012 at 11:51:56AM -0600, Warner Losh wrote: On Jul 8, 2012, at 6:40 AM, David Schultz wrote: On Tue, May 29, 2012, Peter Jeremy wrote: On 2012-May-28 15:54:06 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: Given that cephes was written years before C99 was even conceived, I suspect all functions are sub-standard. Well, most of cephes was written before C99. The C99 parts of cephes were written to turn it into a complete C99 implementation. I'm a bit late to the party, but I thought I'd chime in with some context. We did consider using Cephes years ago, and even got permission from the author to release it under an acceptable license. We later decided not to use it for technical reasons. By the way, virtually none of the people who have complained about the missing functions actually need them. Mostly they just want to compile some software that was written by a naive programmer who thought it would be cool to use the most precise type available. The complex functions are even less commonly needed, and the truth is that they have no business being part of the C standard anyway. The question remains of what to do about the missing functions. Bruce and Steve have been working on expl and logl for years. If those ever get in the tree, the remaining long double functions are easy. Those functions are basically done, modulo a bunch of cleanup and testing, and I encourage any mathematically inclined folks who are interested in pushing things along to get in touch with them. I'm not going to have any time myself for a few months at least. Where can I find these? I've posted expl() a few times for the ld80 version. I don't have an ld128 version, which is why I have yet to submit a formal patch for expl(). I also have an ld80 expm1l(). I have a copy of bde's ld80 logl(). IIRC, bde wrote an ld128, but I don't have nor do I know if it has been tested. Lastly, there's the question of mediocre alternatives, such as solutions that get the boundary cases wrong or don't handle 128-bit floating point. For the exponential and logarithmic functions, Bruce and Steve have already written good implementations, so there's no reason to settle for less. As for the other long double functions, bringing in some Cephes code in a separate directory as a temporary fix might be the way to go. I don't like that solution, and Steve raises some good technical points about why it isn't ideal; however, a better solution is more than a decade overdue, and people are justified in finding that unacceptable. Don't let the perfect be the enemy of the good. It is better to have OK implementations of these functions than none at all. You'll need to define OK, because some of the implementations I've read are poor. I also hope that before anything is committed to libm that the person doing the committing does some rather thorough testing of the code. We originally had so-so double support, but bruce spent many years optimizing them to make them very good. If we were just starting out, and hadn't let 10 years get behind us, I'd give the substandard argument some weight. But now that we're 13 years down the line from c99's publication I think we need to relax our standards a bit. I'd even argue that these functions being a little bad could easily spur people to make them better. Their absence makes people just #define llexp(x) lexp(x), etc. :( Of course, the converse argument is that people have had 10 years to learn enough about floating point to actually write and submit code. I knew very little about FP before I wrote sqrtl(). I had a need for sqrtl() and so I went looking for a solution. PS: I also wrote sincos[fl](), which is very handy for the complex trig functions. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Sun, Jul 08, 2012, Steve Kargl wrote: The question remains of what to do about the missing functions. Bruce and Steve have been working on expl and logl for years. If those ever get in the tree, the remaining long double functions are easy. Those functions are basically done, modulo a bunch of cleanup and testing, and I encourage any mathematically inclined folks who are interested in pushing things along to get in touch with them. I'm not going to have any time myself for a few months at least. Where can I find these? I've posted expl() a few times for the ld80 version. I don't have an ld128 version, which is why I have yet to submit a formal patch for expl(). I also have an ld80 expm1l(). I have a copy of bde's ld80 logl(). IIRC, bde wrote an ld128, but I don't have nor do I know if it has been tested. Yes, Bruce has ld128 versions, and clusteradm very kindly got us a sparc64 machine to test on. That was about the time I ran out of time to keep working on it. If someone wants to pick it up, that would be great. PS: I also wrote sincos[fl](), which is very handy for the complex trig functions. Yes, you should commit that! ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Sun, Jul 08, 2012 at 02:06:46PM -0500, Stephen Montgomery-Smith wrote: Here is a technical question. I know that people always talk about ulp's in the context of how good a function implementation is. I think the ulp is the number of base 2 digits at the end of the mantissa that we cannot rely on. There are a few definitions of ULP (units in the last place). Google 'Muller ULP' and check goldberg's paper (google 'goldberg floating point'. Here's how I define ULP for my testing and the MPFR code. /* From a das email: ulps = err * (2**(LDBL_MANT_DIG - e)), where e = ilogb(z). I use: mpfr_sub(err, acc_out, approx_out, GMP_RNDN); mpfr_abs(err, err, GMP_RNDN); mpfr_mul_2si(ulpx, err, -mpfr_get_exp(acc_out) + LDBL_MANT_DIG, GMP_RNDN); ulp = mpfr_get_d(ulpx, GMP_RNDN); if (ulp 0.506 mpfr_cmpabs(err, ldbl_minx) 0) { print warning...; } */ #include gmp.h #include mpfr.h #include sgk.h static mpfr_t mp_ulp_err; /* * Set the precision of the ULP computation. */ void mp_ulp_init(int prec) { mpfr_init2(mp_ulp_err, (mpfr_prec_t)prec); } /* * Given the 'approximate' value of a function and * an 'accurate' value compute the ULP. */ double mp_ulp(mpfr_t approximate, mpfr_t accurate, int dig) { double ulp; mpfr_sub(mp_ulp_err, accurate, approximate, RD); mpfr_abs(mp_ulp_err, mp_ulp_err, RD); mpfr_mul_2si(mp_ulp_err, mp_ulp_err, - mpfr_get_exp(accurate) + dig, RD); ulp = mpfr_get_d(mp_ulp_err, RD); return (ulp); } So if one were to write a naive implementation of lexp(x) that used Taylor's series if x is positive, and 1/lexp(-x) is x is negative - well one could fairly easily estimate an upper bound on the ulp, but it wouldn't be low (like ulp=1 or 2), but probably rather higher (ulp of the order of 10 or 20). I've already written an ld80 expl(). I have tested billions if not trillions of values. Don't recall the max ULP I saw, but I know it was less than 1. I don't have an ld128 version, so I won't submit a patch for libm. So do people really work hard to get that last drop of ulp out of their calculations? I know very few scientist who work hard to reduce the ULP. Most have little understanding of floating point. Would a ulp=10 be considered unacceptable? Yes, it is unacceptable for the math library. Remember ULPs propagate and can quickly grow if the initial ulp for a result is large. Also, looking through the source code for the FreeBSD implementation of exp, I saw that they used some rather smart rational function. (I don't know how they came up with it.) See Steven Moshier's book, which is the basis for CEPHES. Instead of using a minimax polynomial approximation for the function in some interval, one uses a minimax approximation with a rational function. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 07/08/2012 06:58 PM, Steve Kargl wrote: On Sun, Jul 08, 2012 at 02:06:46PM -0500, Stephen Montgomery-Smith wrote: So do people really work hard to get that last drop of ulp out of their calculations? I know very few scientist who work hard to reduce the ULP. Most have little understanding of floating point. Would a ulp=10 be considered unacceptable? Yes, it is unacceptable for the math library. Remember ULPs propagate and can quickly grow if the initial ulp for a result is large. OK. But suppose I wanted ld80 precision. What would be the advantage of using an ld80 expl function with a ulp of 1 over an ld128 expl function with a ulp of 10? The error propagation in the latter case could not be worse than the error propagation in the latter case. In other words, if I were asked to write a super-fantastic expl function, where run time was not a problem, I would use mpfr, use Taylor's series with a floating point precision that had way more digits than I needed, and then just truncate away the last digits when returning the answer. And this would be guaranteed to give the correct answer to just above 0.5 ulp (which I assume is best possible). From a scientist's point of view, I would think ulp is a rather unimportant concept. The concepts of absolute and relative error are much more important to them. The only way I can see ULP errors greatly propagating is if one is performing iteration type calculations (like f(f(f(f(x). This sort of thing is done when computing Julia sets and such like. And in that case, as far as I can see, a slightly better ulp is not going to drastically change the limitations of whatever floating point precision you are using. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Sun, Jul 08, 2012 at 07:29:30PM -0500, Stephen Montgomery-Smith wrote: On 07/08/2012 06:58 PM, Steve Kargl wrote: On Sun, Jul 08, 2012 at 02:06:46PM -0500, Stephen Montgomery-Smith wrote: So do people really work hard to get that last drop of ulp out of their calculations? I know very few scientist who work hard to reduce the ULP. Most have little understanding of floating point. Would a ulp=10 be considered unacceptable? Yes, it is unacceptable for the math library. Remember ULPs propagate and can quickly grow if the initial ulp for a result is large. OK. But suppose I wanted ld80 precision. What would be the advantage of using an ld80 expl function with a ulp of 1 over an ld128 expl function with a ulp of 10? The error propagation in the latter case could not be worse than the error propagation in the latter case. Well, on the most popular hardware (that being i386/amd64), ld80 will use hardware fp instruction while ld128 must be done completely in software. The speed difference is significant. In other words, if I were asked to write a super-fantastic expl function, where run time was not a problem, I would use mpfr, use Taylor's series with a floating point precision that had way more digits than I needed, and then just truncate away the last digits when returning the answer. And this would be guaranteed to give the correct answer to just above 0.5 ulp (which I assume is best possible). It's more like 1 ULP after truncation (, which isn't truncation but rounding). The problem is run-time that 'run-time is the problem'. Try writing a expl() and simply use mpfr_exp() with 64-bit precision. If you're doing any serious simulation where exp() will be evaluate millions or billions of time, you'll notice the difference. The only way I can see ULP errors greatly propagating is if one is performing iteration type calculations (like f(f(f(f(x). Have you read Goldberg's paper? Not to mention, I've seen way too many examples of 'x - y' where cancellation of significant digits causes problems. Throw in rather poor estimates of function results with real poor ULP and you have problems. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Jul 8, 2012, at 8:01 PM, Steve Kargl wrote: Not to mention, I've seen way too many examples of 'x - y' where cancellation of significant digits causes problems. Throw in rather poor estimates of function results with real poor ULP and you have problems. Are these problems significantly more or less than the usual #define I talked about before? If the functions are so so, but much better than the double version, we have a significant win, even if things aren't perfect. If we weren't 13 past the publication date of the c99 standard, I'd be more sympathetic to the 'we need a high quality implementation' arguments. However, we can't let the perfect be the enemy of the good here. We claim c99 conformance, yet don't have these functions. After all, many of the original functions that were in our library had sub-optimial performance which bruce optimized over many years. Why can't we use this model here? Warner ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Sun, Jul 08, 2012 at 08:13:21PM -0600, Warner Losh wrote: On Jul 8, 2012, at 8:01 PM, Steve Kargl wrote: Not to mention, I've seen way too many examples of 'x - y' where cancellation of significant digits causes problems. Throw in rather poor estimates of function results with real poor ULP and you have problems. Are these problems significantly more or less than the usual #define I talked about before? If the functions are so so, but much better than the double version, we have a significant win, even if things aren't perfect. The answer is more complicated than the question asked. On i386, the precision of long double and double is 53 bits. On amd64, the precision of long double and double are 64 and 53 bits. On i386 and amd64, emax and emin are 1024 and -1023 (could be off by one here) for double and 16384 and -16383 for long double. Now, consider 2**emin = exp(x) = 2**emax for some set of x (where I'm using the Fortran exponential operator ** instead of powl(2,{emin|emax}) because FreeBSD does not have a powl()). The sets for long double and double are significantly different. If we weren't 13 past the publication date of the c99 standard, I'd be more sympathetic to the 'we need a high quality implementation' arguments. However, we can't let the perfect be the enemy of the good here. We claim c99 conformance, yet don't have these functions. AFAIK, neither gcc in base nor clang would be c99 complaint even if all of the c99 math functions were available. After all, many of the original functions that were in our library had sub-optimial performance which bruce optimized over many years. Why can't we use this model here? I think the fear is that committing sub-optimal code would take many many more years to fix. AFAICT, bde longer commits code to src/ (since the switch from cvs to svn), and he depends on his code reviews to get things fixed. With das and me both having little to no time for FreeBSD, that leaves no one to review/commit his suggested fixes (not that I would question bde's opinion on a code change). -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 07/08/2012 09:01 PM, Steve Kargl wrote: Have you read Goldberg's paper? I must admit that I had not. I found it at: http://docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html Not to mention, I've seen way too many examples of 'x - y' where cancellation of significant digits causes problems. Throw in rather poor estimates of function results with real poor ULP and you have problems. I think that if a program includes x-y where cancellation causes huge loss of digits, then the program should be considered highly flawed. The programmer might get lucky, and a few extra ULP might save his skin, but it would be better if the program catastrophically failed so that he would realize they need to do something smarter. I got my first scientific calculator around 1975, and I remember my surprise when I discovered that (1/3)*3-1 on my calculator did not produce zero. Years later I bought another calculator, and imagine my further surprise when (1/3)*3-1 did produce zero. They must have done something very smart, I thought. The problem is, these kinds of tricks can only help you so much. Calculations like arccos(arccos(arccos(arccos(cos(cos(cos(cos(x)-x are probably not going to be correct no matter how smart the programmer is. The paper by Goldberg has some really nice stuff. I teach a class on numerical methods. One example I present is the problem using equation (4) of Goldberg's paper to solve quadratic equations. I tell them that the smart thing to do is to use equation (4) or equation (5) depending on whether b has the same sign as sqrt(b^2-4ac). This is a very good illustration of how to overcome the x-y problem, and in my opinion it has to be understood by the programmer, not the writer of the square-root package. Trying to compensate by getting that lost drop of ULP out of square root is looking in the wrong direction. But there is other stuff in his paper I don't like, and it comes across as nit-picking rather than really thinking through why he wants the extra accuracy. I feel like he is saving the pennies, but losing the dollars. For example, computing the quadratic formula when b^2 and 4ac are roughly equal (discussed in his proof of Theorem 4). He says that roughly half the digits of the final answer will be contaminated. He recommends performing the calculation with double the precision, and then retaining what is left. The reason I don't like his solution is this. The contamination of half the digits is real, not some artifact of the computing method. Unless you know EXACTLY what a, b and c are, only half the digits of the quadratic formula are valid, no matter how you calculate it. For example, maybe a is calculated by some other part of the program as 1/3. Since 1/3 cannot be represented exactly in base 2, it will not be exactly 1/3. That part of the program doesn't know that this number is going to be used later on in the quadratic formula, so it computes it in single precision. Now the quadratic formula algorithm converts a to double precision. But a will still be 1/3 to single precision. The only way to avoid this error is to perform EVERY pertinent calculation prior to the use of quadratic formula in double precision. (And suppose instead the data comes from a scientific experiment - the programmer needs to be performing an error analysis, not hoping his implementation of quadratic formula is performed in double precision.) Similarly, I am not a fan of the Kahan Summation Formula (Theorem 8). There are two reasons why I might compute large sums. One is accounting, when I better be using high enough precision that the arithmetic is exact. The other is something like numerical integration. But in the latter case, unless the summands have a very special form, it is likely that each summand will only be accurate to within e. Thus the true sum is only known to within n*e, and so the naive method really hasn't lost any precision at all. And typically n*e will be so small compared to the true answer that it won't matter. If it does matter, the problem is that n is too big. The algorithm will take decades to run. Focus efforts on producing a different numerical integration method, instead of getting that lost drop of ULP. I cannot envision any non-contrived circumstance where the Kahan summation formula is going to give an answer that is significantly more reliable than naive summation. I can envision a circumstance when the reduction in speed of the Kahan summation formula over naive summation could be significant. I think the example I like best that illustrates floating point errors is inverting badly conditioned matrices. And any scientist who works with large matrices had better know what ill-conditioned means. The proper answer is that if your matrix is ill-conditioned, give up. You need to look at the problem where your matrix came from, and rethink how you concert
Re: Use of C99 extra long double math functions after r236148
On Sun, Jul 08, 2012 at 10:41:44PM -0500, Stephen Montgomery-Smith wrote: On 07/08/2012 09:01 PM, Steve Kargl wrote: Have you read Goldberg's paper? I must admit that I had not. I found it at: http://docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html Yes, it's easy to find. Unfortunately, many people don't know that it exists. I've read that paper between 10 and 20 times and I still learn something new with every reading. Not to mention, I've seen way too many examples of 'x - y' where cancellation of significant digits causes problems. Throw in rather poor estimates of function results with real poor ULP and you have problems. I think that if a program includes x-y where cancellation causes huge loss of digits, then the program should be considered highly flawed. The programmer might get lucky, and a few extra ULP might save his skin, but it would be better if the program catastrophically failed so that he would realize they need to do something smarter. I got my first scientific calculator around 1975, and I remember my surprise when I discovered that (1/3)*3-1 on my calculator did not produce zero. Years later I bought another calculator, and imagine my further surprise when (1/3)*3-1 did produce zero. They must have done something very smart, I thought. Yes. They used CORDIC arithmetic instead of (binary) floating point. The problem is, these kinds of tricks can only help you so much. Calculations like arccos(arccos(arccos(arccos(cos(cos(cos(cos(x)-x are probably not going to be correct no matter how smart the programmer is. Well, I would claim that any programmer who wrote such an expression is not smart. The paper by Goldberg has some really nice stuff. I teach a class on numerical methods. It may be profitable to at least have your students read the paper. I don't know about others, but I find doing floating point arithematic correctly very difficult. One example I present is the problem using equation (4) of Goldberg's paper to solve quadratic equations. I tell them that the smart thing to do is to use equation (4) or equation (5) depending on whether b has the same sign as sqrt(b^2-4ac). This is a very good illustration of how to overcome the x-y problem, and in my opinion it has to be understood by the programmer, not the writer of the square-root package. Trying to compensate by getting that lost drop of ULP out of square root is looking in the wrong direction. Careful. IEEE 754 specifies and one can prove that sqrt() can be computed correctly to less than or equal to 1/2 ULP for all input in all rounding modes. Computing the roots of the quartic equation, which uses the sqrt() function, is a different matter. But there is other stuff in his paper I don't like, and it comes across as nit-picking rather than really thinking through why he wants the extra accuracy. I feel like he is saving the pennies, but losing the dollars. well, in Goldberg's defense, think about the title of the paper. I think his major aim in the paper is get scientist to think about what and how they compute some number. Similarly, I am not a fan of the Kahan Summation Formula (Theorem 8). Oh my. I'm a fan of Kahan's summation formula. But, I assume that one applies it to situations where it is needed. There are two reasons why I might compute large sums. One is accounting, when I better be using high enough precision that the arithmetic is exact. If you're doing accouting, hopefully, you're using BCD. Accouting and floating point simply do not mix. The other is something like numerical integration. But in the latter case, unless the summands have a very special form, it is likely that each summand will only be accurate to within e. Thus the true sum is only known to within n*e, and so the naive method really hasn't lost any precision at all. And typically n*e will be so small compared to the true answer that it won't matter. If it does matter, the problem is that n is too big. The algorithm will take decades to run. Focus efforts on producing a different numerical integration method, instead of getting that lost drop of ULP. I cannot envision any non-contrived circumstance where the Kahan summation formula is going to give an answer that is significantly more reliable than naive summation. I can envision a circumstance when the reduction in speed of the Kahan summation formula over naive summation could be significant. I can envision a use. Try computing the RMS height of a numerically generated random rough surface in the scattering of an incident acoustic field from said surface. The RMS height is a main component in the first order perturbation theory expansion of the scattered field. You get the RMS height wrong, you have no idea what the scatter field strength is. I think the example I like best that illustrates floating point errors is inverting badly
Re: Use of C99 extra long double math functions after r236148
On 31 May 2012 08:45, John Baldwin j...@freebsd.org wrote: I do think we should provide something in ports as an interim solution. There are other 3rd party applications looking to drop FreeBSD support because we are missing APIs that almost all other OS's have. Â I'm fine if the interim lives in ports and that we don't import substandard routines into the base. Â I would even be fine with calling it /usr/local/lib/libm_inaccurate.so. Â However, I do think we need an option. Do we have a wiki page listing the functions in libm we are missing? Having some kind of place to track progress and figure out what exactly is needed is the first step to getting these APIs into shape. Also, are there BSD licensed naive implementations of these functions we can use? Would it be okay to have slow, but accurate versions of these functions as a stopgap? -- Eitan Adler ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
I do think we should provide something in ports as an interim solution. There are other 3rd party applications looking to drop FreeBSD support because we are missing APIs that almost all other OS's have. I'm fine if the interim lives in ports and that we don't import substandard routines into the base. I would even be fine with calling it /usr/local/lib/libm_inaccurate.so. However, I do think we need an option. I think it should be called libm.so. Otherwise we have to do a serious editing job on the Makefiles/configure scripts. I think that this is as it should be: only those applications that really need to link against such a library should do so, and this should be enforced and checked by the port maintainer. If you call it by the same name you still need to make other changes to ensure that the right library is used, and these other changes can be confusing and lead to problems, as with the openssl and gcc support libraries. Some portable applications already provide convenient variables in their build infrastructure, since the quality, coverage, and location of the math libraries varies on different systems. b. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
Do we have a wiki page listing the functions in libm we are missing? Having some kind of place to track progress and figure out what exactly is needed is the first step to getting these APIs into shape. I already suggested this, and mentioned: http://wiki.freebsd.org/MissingMathStuff Also, are there BSD licensed naive implementations of these functions we can use? Would it be okay to have slow, but accurate versions of these functions as a stopgap? I don't know of any BSD-licensed routines that would be suitable for use in the base system. You can cobble together replacements from various ports, with different licenses. For example, if you wanted accuracy and correctness, but didn't care about speed, you could use math/mpc and math/mpfr. Unfortunately, with many naive implementations, you often lose more than just speed. b. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Friday, June 01, 2012 1:55:10 am Eitan Adler wrote: On 31 May 2012 08:45, John Baldwin j...@freebsd.org wrote: I do think we should provide something in ports as an interim solution. There are other 3rd party applications looking to drop FreeBSD support because we are missing APIs that almost all other OS's have. I'm fine if the interim lives in ports and that we don't import substandard routines into the base. I would even be fine with calling it /usr/local/lib/libm_inaccurate.so. However, I do think we need an option. Do we have a wiki page listing the functions in libm we are missing? Having some kind of place to track progress and figure out what exactly is needed is the first step to getting these APIs into shape. Also, are there BSD licensed naive implementations of these functions we can use? Would it be okay to have slow, but accurate versions of these functions as a stopgap? Peter Jeremy more or less has a stopgap already ready judging by the comments in the thread thus far. -- John Baldwin ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 2012-Jun-01 10:29:13 -0400, John Baldwin j...@freebsd.org wrote: On Friday, June 01, 2012 1:55:10 am Eitan Adler wrote: Also, are there BSD licensed naive implementations of these functions we can use? Would it be okay to has slow, but accurate versions of these functions as a stopgap? Peter Jeremy more or less has a stopgap already ready judging by the comments in the thread thus far. There's probably an hours work by either stephen@ or myself to adapt the work I did on cephes in Sage to a standalone FreeBSD port. Unfortunately, both stephen@ I are currently otherwise occupied and other comments in this thread suggest that the inclusion of such a port would be strongly opposed. Note that cephes isn't slow but accurate - it's reasonably fast but naive and therefore dodgy in edge cases. -- Peter Jeremy pgpHAsPC0mWbI.pgp Description: PGP signature
Re: Use of C99 extra long double math functions after r236148
On 1 June 2012 17:03, Peter Jeremy pe...@rulingia.com wrote: On 2012-Jun-01 10:29:13 -0400, John Baldwin j...@freebsd.org wrote: On Friday, June 01, 2012 1:55:10 am Eitan Adler wrote: Also, are there BSD licensed naive implementations of these functions we can use? Would it be okay to has slow, but accurate versions of these functions as a stopgap? Peter Jeremy more or less has a stopgap already ready judging by the comments in the thread thus far. There's probably an hours work by either stephen@ or myself to adapt the work I did on cephes in Sage to a standalone FreeBSD port. Unfortunately, both stephen@ I are currently otherwise occupied and other comments in this thread suggest that the inclusion of such a port would be strongly opposed. Note that cephes isn't slow but accurate - it's reasonably fast but naive and therefore dodgy in edge cases. Yes, I was asking if any of the former type exist. Optimally we would want fast and accurate - but it doesn't currently exist. Fast, but inaccurate has been strongly objected to. Is there third option? -- Eitan Adler ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Fri, Jun 01, 2012 at 05:16:03PM -0700, Eitan Adler wrote: On 1 June 2012 17:03, Peter Jeremy pe...@rulingia.com wrote: On 2012-Jun-01 10:29:13 -0400, John Baldwin j...@freebsd.org wrote: On Friday, June 01, 2012 1:55:10 am Eitan Adler wrote: Also, are there BSD licensed naive implementations of these functions we can use? Would it be okay to has slow, but accurate versions of these functions as a stopgap? Peter Jeremy more or less has a stopgap already ready judging by the comments in the thread thus far. There's probably an hours work by either stephen@ or myself to adapt the work I did on cephes in Sage to a standalone FreeBSD port. Unfortunately, both stephen@ I are currently otherwise occupied and other comments in this thread suggest that the inclusion of such a port would be strongly opposed. Note that cephes isn't slow but accurate - it's reasonably fast but naive and therefore dodgy in edge cases. Yes, I was asking if any of the former type exist. Optimally we would want fast and accurate - but it doesn't currently exist. Fast, but inaccurate has been strongly objected to. Is there third option? Of course. Sit down and write code. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 1 June 2012 23:52, Steve Kargl s...@troutmask.apl.washington.edu wrote: Of course. Â Sit down and write code. If I ever find the time, I just might. Do we have a wiki page listing the set of functions which we don't yet have? -- Eitan Adler ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Sat, Jun 02, 2012 at 12:04:58AM -0400, Eitan Adler wrote: On 1 June 2012 23:52, Steve Kargl s...@troutmask.apl.washington.edu wrote: Of course. ??Sit down and write code. If I ever find the time, I just might. Do we have a wiki page listing the set of functions which we don't yet have? I don't know. I don't read the wiki. You can many/most/all by looking in /usr/include/math.h There are many #if 0 ... #endif blocks. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Monday, May 28, 2012 7:02:18 pm Steve Kargl wrote: On Tue, May 29, 2012 at 07:05:07AM +1000, Peter Jeremy wrote: On 2012-May-28 11:01:24 -0500, Stephen Montgomery-Smith step...@missouri.edu wrote: One thing that could be done is to have a math/cephes port that adds the extra C99 math functions. This is already done in the math/sage port, using a rather clever patch due to Peter Jeremy, that applies to the cephes code. What it would do is to create a /usr/local/lib/libm.so that would provide the extra functions not currently included in /lib/libm.so, and then link in /lib/libm.so as well. It would also create its own /usr/local/include/math.h and /usr/local/include/complex.h as well. Basically, as long as the compiler searches /usr/local/{include,lib} before the base include/lib then math.h, complex.h and -lm give the application a complete C99 math implementation by using base functions where they exist and cephes functions where they don't. The patch I wrote for sage can be found at http://trac.sagemath.org/sage_trac/ticket/9543 If there's any interest, I could produce a port for this. Another option would be to import cephes into base and use it to provide the missing C99 functions. Cephes includes copyright notices but the closest I can find to a license is: Some software in this archive may be from the book _Methods and Programs for Mathematical Functions_ (Prentice-Hall or Simon Schuster International, 1989) or from the Cephes Mathematical Library, a commercial product. In either event, it is copyrighted by the author. What you see here may be used freely but it comes with no support or guarantee. Please talk to das@ (although I believe he's finishing up his dissertation). I recall that he's stated that he looked into using cephes, and concluded that it is not suitable for libm. Note there is also http://www.freebsd.org/cgi/query-pr.cgi?pr=kern/147599 which I've also objected to importing into libm. I do think we should provide something in ports as an interim solution. There are other 3rd party applications looking to drop FreeBSD support because we are missing APIs that almost all other OS's have. I'm fine if the interim lives in ports and that we don't import substandard routines into the base. I would even be fine with calling it /usr/local/lib/libm_inaccurate.so. However, I do think we need an option. -- John Baldwin ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 05/31/12 10:45, John Baldwin wrote: On Monday, May 28, 2012 7:02:18 pm Steve Kargl wrote: On Tue, May 29, 2012 at 07:05:07AM +1000, Peter Jeremy wrote: On 2012-May-28 11:01:24 -0500, Stephen Montgomery-Smith step...@missouri.edu wrote: One thing that could be done is to have a math/cephes port that adds the extra C99 math functions. This is already done in the math/sage port, using a rather clever patch due to Peter Jeremy, that applies to the cephes code. What it would do is to create a /usr/local/lib/libm.so that would provide the extra functions not currently included in /lib/libm.so, and then link in /lib/libm.so as well. It would also create its own /usr/local/include/math.h and /usr/local/include/complex.h as well. Basically, as long as the compiler searches /usr/local/{include,lib} before the base include/lib thenmath.h,complex.h and -lm give the application a complete C99 math implementation by using base functions where they exist and cephes functions where they don't. The patch I wrote for sage can be found at http://trac.sagemath.org/sage_trac/ticket/9543 If there's any interest, I could produce a port for this. Another option would be to import cephes into base and use it to provide the missing C99 functions. Cephes includes copyright notices but the closest I can find to a license is: Some software in this archive may be from the book _Methods and Programs for Mathematical Functions_ (Prentice-Hall or Simon Schuster International, 1989) or from the Cephes Mathematical Library, a commercial product. In either event, it is copyrighted by the author. What you see here may be used freely but it comes with no support or guarantee. Please talk to das@ (although I believe he's finishing up his dissertation). I recall that he's stated that he looked into using cephes, and concluded that it is not suitable for libm. Note there is also http://www.freebsd.org/cgi/query-pr.cgi?pr=kern/147599 which I've also objected to importing into libm. I do think we should provide something in ports as an interim solution. There are other 3rd party applications looking to drop FreeBSD support because we are missing APIs that almost all other OS's have. I'm fine if the interim lives in ports and that we don't import substandard routines into the base. I would even be fine with calling it /usr/local/lib/libm_inaccurate.so. However, I do think we need an option. I think it should be called libm.so. Otherwise we have to do a serious editing job on the Makefiles/configure scripts. sed -E 's/[[::]]-lm[[::]]/-lm_inaccuarate/' might have some false positives and false negatives. (Did I even get the sed syntax correct?) By the way, I was looking through: http://www.freebsd.org/cgi/query-pr.cgi?pr=kern/147599 They say they compute cosh(2y) - cos(2x) using Taylor's series. I would agree with Steve that this is seriously low quality. Floating point errors are likely to be truly horrible. This will fail even for non-edge cases. (This was the same error that openoffice used to have for computing erf(x), that gave things like erf(100) = -4353243242.) ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 01/06/2012, at 2:42, Stephen Montgomery-Smith wrote: I do think we should provide something in ports as an interim solution. There are other 3rd party applications looking to drop FreeBSD support because we are missing APIs that almost all other OS's have. I'm fine if the interim lives in ports and that we don't import substandard routines into the base. I would even be fine with calling it /usr/local/lib/libm_inaccurate.so. However, I do think we need an option. I think it should be called libm.so. Otherwise we have to do a serious editing job on the Makefiles/configure scripts. sed -E 's/[[::]]-lm[[::]]/-lm_inaccuarate/' might have some false positives and false negatives. (Did I even get the sed syntax correct?) Another option would be to put it in base but bleat about it if it's actually used (like mktemp et al) -- Daniel O'Connor software and network engineer for Genesis Software - http://www.gsoft.com.au The nice thing about standards is that there are so many of them to choose from. -- Andrew Tanenbaum GPG Fingerprint - 5596 B766 97C0 0E94 4347 295E E593 DC20 7B3F CE8C ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 05/29/12 19:54, Stephen Montgomery-Smith wrote: [...] Anyway, given that floating point is a big issue, and we are about a decade behind schedule, really suggests that a floating-po...@freebsd.org mailing list is needed. Or maybe there is an existing freebsd mailing list you guys already occupy. I would certainly take part in discussions on such a list. And since I have developed similar software in the past, I could even contribute a bit more, if desired. I have quite a tight schedule, but this type of development does not require an extensive infrastructure and can be done in lost hours on the road. Kind regards, Hans Ottevanger ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
Perhaps a more general name might be appropriate so as to include fixed-point problems? math-libs@? numerics@? I'm sure someone will come up with a better name. mcl ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Wed, May 30, 2012 at 8:51 AM, Mark Linimon lini...@lonesome.com wrote: Perhaps a more general name might be appropriate so as to include fixed-point problems? math-libs@? numerics@? I'm sure someone will come up with a better name. mcl Numerics@ is good . Other names may be numerical-analysis@ computing@ ( includes high performance computing ) mathematics@ ( includes statistics , numerical analysis , computation , probability , ... ) Thank you very much . Mehmet Erol Sanliturk ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
Holy crap, if there was ever a current example of a bikeshed, this is it. :-) Adrian On 30 May 2012 10:19, Mehmet Erol Sanliturk m.e.sanlit...@gmail.com wrote: On Wed, May 30, 2012 at 8:51 AM, Mark Linimon lini...@lonesome.com wrote: Perhaps a more general name might be appropriate so as to include fixed-point problems? Â math-libs@? Â numerics@? I'm sure someone will come up with a better name. mcl Numerics@ is good . Other names may be numerical-analysis@ computing@ Â Â Â Â Â Â Â ( includes high performance computing ) mathematics@ Â Â Â Â Â ( includes statistics , numerical analysis , computation , probability , ... ) Thank you very much . Mehmet Erol Sanliturk ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
This discussion confirms my impression, that it should be possible as an interim solution, to use a port for missing math functions (cephes alike or whatever). The port itself could warn the user about inaccuracies and edge-cases. Parts of Cephes are already in ports: math/ldouble. I had planned to add the remainder. It can be useful, but it has problems, some of which have been described. The same is true of other open-source alternatives that I am aware of. I would be happy to take on any steep learning curve, and make contributions, if only I were part of a group that would steer me in the right direction. A Functional Analyst offering to write code for a floating-point math library!?!! ;) You should send me a message off-list: maybe I can find some time to help. Anyway, given that floating point is a big issue, and we are about a decade behind schedule, really suggests that a floating-point at freebsd.org mailing list is needed. Or maybe there is an existing freebsd mailing list you guys already occupy. I do not know that a separate FreeBSD mailing list would be of great help, considering the likely volume of traffic, and the fact that we already have the standards mailing list. There is also: http://mailman.oakapple.net/pipermail/numeric-interest/ which serves, among other things, as a forum for discussing changes and improvements to the open-source math library from which parts of our system math library were derived: http://mailman.oakapple.net/pipermail/numeric-interest/2010-September/002054.html A wiki page detailing problems and procedures, with a wish list for the system libraries and toolchain(s) might help. At present there are only scattered messages in the mailing list archives, PRs, and http://wiki.freebsd.org/MissingMathStuff . b. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Tue, May 29, 2012 at 02:56:13PM +1000, Peter Jeremy wrote: On 2012-May-28 15:54:06 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: There some test code in cephes. Can you point me to a suitable test suite for LD80 and LD128? The reason for calling it libm is to avoid having to hack every consumer to add an additional library. I can't point you to test code. When I work on a function, I write test code. It took me 3+ years to get sqrtl() into libm, but bde and das (and myself) wanted to make sure the code worked. Last time I checked (a couple of years ago), FreeBSD was missing 65 C99 libm functions. At 3 years per function, we should have C99 support available early in the 23rd century - which may be a bit late. sqrtl() is a bit special in that IEEE 754 requires that it have no more than 0.5 ULP for all arguments in all roundng modes. As to other functions, I've been trying for 10+ years to get some of these into FreeBSD. I can assure you that there has never been a rush of people volunteering to help or a rush of people willing to fund the development of the necessary code. On 2012-May-28 22:03:43 -0500, Stephen Montgomery-Smith step...@missouri.edu wrote: 1. By being so picky about being so precise, FreeBSD is behind the time line in rolling out a usable set of C99 functions. And at the current rate, we'll all be long dead before they are available. It seems you've had a couple of years to implement one or more of the missing functions. Yes, we'll all be dead before libm is C99 ready because no one, other than bde@, das@ and myself, has been willing to actually help write the code. I know that at least one of those three people has had no time in the last year or two work on libm. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 29.05.2012 08:10 (UTC+1), Steve Kargl wrote: On Tue, May 29, 2012 at 02:56:13PM +1000, Peter Jeremy wrote: On 2012-May-28 15:54:06 -0700, Steve Kargls...@troutmask.apl.washington.edu wrote: There some test code in cephes. Can you point me to a suitable test suite for LD80 and LD128? The reason for calling it libm is to avoid having to hack every consumer to add an additional library. I can't point you to test code. When I work on a function, I write test code. It took me 3+ years to get sqrtl() into libm, but bde and das (and myself) wanted to make sure the code worked. Last time I checked (a couple of years ago), FreeBSD was missing 65 C99 libm functions. At 3 years per function, we should have C99 support available early in the 23rd century - which may be a bit late. sqrtl() is a bit special in that IEEE 754 requires that it have no more than 0.5 ULP for all arguments in all roundng modes. As to other functions, I've been trying for 10+ years to get some of these into FreeBSD. I can assure you that there has never been a rush of people volunteering to help or a rush of people willing to fund the development of the necessary code. The almost complete absence of volunteers in this case is first of all caused by the need of very special, deep knowlegdes on mathematical and informatical issues. I bet there are only a few people out there, who are really good in both (three of them are known to us ;) ). The problems with quality standards in such libraries spreaded over most OS'es and platforms support this view. (I am aware that these arguments are kown and discussed before.) People like me, who wanted to _use_ FreeBSD as scientific desktops, have an urgent need on such mathematical libraries. But most of us have no chance to contribute conceptual. For some tasks we could use specialized libraries or software, which include missing procedures and functions. But more often, modern apps expect to get this functionalities from the OS. This discussion confirms my impression, that it should be possible as an interim solution, to use a port for missing math functions (cephes alike or whatever). The port itself could warn the user about inaccuracies and edge-cases. The more important question to me is, how can the remaining (huge!) work on the systems libm done in a foreseeable time. As one possibility, wouldn't it be imaginable to pay some people for doing (some of) this work (FreeBSD Foundation, Sponsors from industry and science etc.)? I personally, as a private person, would be willing to pay a few hundred dollars for it. Paying people for other tasks in FreeBSD is not entirely uncommon, if I am not mistaken. On 2012-May-28 22:03:43 -0500, Stephen Montgomery-Smithstep...@missouri.edu wrote: 1. By being so picky about being so precise, FreeBSD is behind the time line in rolling out a usable set of C99 functions. And at the current rate, we'll all be long dead before they are available. It seems you've had a couple of years to implement one or more of the missing functions. Yes, we'll all be dead before libm is C99 ready because no one, other than bde@, das@ and myself, has been willing to actually help write the code. I know that at least one of those three people has had no time in the last year or two work on libm. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 05/29/2012 11:48 AM, Rainer Hurling wrote: On 29.05.2012 08:10 (UTC+1), Steve Kargl wrote: sqrtl() is a bit special in that IEEE 754 requires that it have no more than 0.5 ULP for all arguments in all roundng modes. As to other functions, I've been trying for 10+ years to get some of these into FreeBSD. I can assure you that there has never been a rush of people volunteering to help or a rush of people willing to fund the development of the necessary code. The almost complete absence of volunteers in this case is first of all caused by the need of very special, deep knowlegdes on mathematical and informatical issues. I bet there are only a few people out there, who are really good in both (three of them are known to us ;) ). The problems with quality standards in such libraries spreaded over most OS'es and platforms support this view. (I am aware that these arguments are kown and discussed before.) The more important question to me is, how can the remaining (huge!) work on the systems libm done in a foreseeable time. As one possibility, wouldn't it be imaginable to pay some people for doing (some of) this work (FreeBSD Foundation, Sponsors from industry and science etc.)? I personally, as a private person, would be willing to pay a few hundred dollars for it. Paying people for other tasks in FreeBSD is not entirely uncommon, if I am not mistaken. If there were a freebsd-floating-po...@freebsd.org mailing list to discuss these issues, I would eagerly join. While I have never written libm-style libraries myself, I do have knowledge of both computers AND mathematics, and I do teach two university courses: Numerical Methods and Numerical Linear Algebra. And I have written quite a bit of numerical code, so I do know about the issues of testing and debugging such code. I would be happy to take on any steep learning curve, and make contributions, if only I were part of a group that would steer me in the right direction. Finally, Steve made a point that the way gcc multiplies complex numbers is wasteful in computational time. I have been bitten by this. I wrote some code that involved multiplying large numbers of complex numbers. I first wrote it using the complex C99 code. Then I rewrote it using double, writing out the complex multiplication long hand, and the program went twice as fast! (It involved multiplying numbers that were purely real or purely imaginary, so performing (x)*(I*y) by (x+0*I)*(0+y*I) really did slow it down.) Anyway, given that floating point is a big issue, and we are about a decade behind schedule, really suggests that a floating-po...@freebsd.org mailing list is needed. Or maybe there is an existing freebsd mailing list you guys already occupy. Stephen ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
.. I don't think it's a bad idea to get freebsd-numeric@ created and start discussing this kind of stuff there. Adrian ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
Hi, Reference: From: Stephen Montgomery-Smith step...@missouri.edu Anyway, given that floating point is a big issue, and we are about a decade behind schedule, really suggests that a floating-po...@freebsd.org mailing list is needed. Or maybe there is an existing freebsd mailing list you guys already occupy. The string fl does not occur in http://lists.freebsd.org/mailman/listinfo/ Apart from this list, you might also find extra people interested to support a proposal to create a floating-point@ list among the subscribers to: freebsd-performance@ freebsd-standards@ freebsd-toolchain@ Good luck Cheers, Julian -- Julian Stacey, BSD Unix Linux C Sys Eng Consultants Munich http://berklix.com Reply below not above, cumulative like a play script, indent with . Format: Plain text. Not HTML, multipart/alternative, base64, quoted-printable. Mail from @yahoo dumped @berklix. http://berklix.org/yahoo/ ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 05/28/2012 07:07 PM, Steve Kargl wrote: On Mon, May 28, 2012 at 06:44:42PM -0500, Stephen Montgomery-Smith wrote: On 05/28/2012 06:30 PM, Steve Kargl wrote: From clog.c in http://www.netlib.org/cephes/c9x-complex double complex ccosh (z) double complex z; { double complex w; double x, y; x = creal(z); y = cimag(z); w = cosh (x) * cos (y) + (sinh (x) * sin (y)) * I; return (w); } See math_private.h about the above. I looked in math_private.h - I presume you meant lib/msun/src/math_private.h. I wasn't able to find anything about ccosh there. I think that for a rough and ready ccosh, this is high enough quality for a math/cephes port. That's the problem. It is not acceptable (imo). The comment in math_private.h that is relevant is /* * Inline functions that can be used to construct complex values. * * The C99 standard intends x+I*y to be used for this, but x+I*y is * currently unusable in general since gcc introduces many overflow, * underflow, sign and efficiency bugs by rewriting I*y as * (0.0+I)*(y+0.0*I) and laboriously computing the full complex product. * In particular, I*Inf is corrupted to NaN+I*Inf, and I*-0 is corrupted * to -0.0+I*0.0. */ Those wrong +-0 mean you may up end up on the worng riemann sheet, and that NaN propagates. OK, I agree with you that gcc fails to be C99 compliant. But I disagree with you that this is a big deal. 1. By being so picky about being so precise, FreeBSD is behind the time line in rolling out a usable set of C99 functions. The compilers and operating systems have many bugs, and if we waited until we were totally sure that all the bugs were gone, we wouldn't have any operating system to work with at all. Why be more picky with floating point arithmetic than the other aspects of FreeBSD? 2. If I was really worried about being on the correct Riemann sheet, I would code very, very carefully, and not rely on numerical accuracy of any kind of floating point calculation. There is a certain mathematical inconsistency in introducing 0, -0 and Inf+I*Inf, and the number of times programmers really want a very specific kind of behavior for exceptions is so rare that they are better off writing their own wrapper code than relying on any library functions. (For example, is there a difference between 0 and +0? If not, does -(0) compute to 0 or -0? I can see circumstances where I sometimes want one, and sometimes the other.) 3. But to counter my own argument, it highly bothers me that in C that (-5)%3 evaluates to -2 and not to 1. That bug^H^H^H feature has truly bitten me. And I have had lengthy arguments online with some C experts as to why 1 should be the correct answer, without being able to convince them. If it were up to me, the whole of the C standard would be scrapped, and rewritten with everything exactly the same except for this one thing. But to those I argued with, I seem just as picky as you seem when you insist that 0 and -0 are different. So what do I do? I live with it, just like we all live in an imperfect world. If I had a choice between correcting the C standard for % or solving world hunger, I would definitely settle for the second. And maybe the C programming language, for all its imperfections, has helped push frontiers of technology sufficiently that someone is less hungry than they would have otherwise been. And if those resources used to feed people had been redirected to fix the C standard, then maybe a few more people would be hungry. In the end, I do think it is good to ultimately settle on good C99 compliant code. But having something intermediate that mostly works is better than nothing. Especially if it exists only in the ports, and not in the base code. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 28 May 2012, at 05:35, Rainer Hurling wrote: Yesterday r236148 (Allow inclusion of libc++ cmath to work after including math.h) was comitted to head, many thanks. Does this mean, that extra long double functions like acoshl, expm1l or log1pl are now really implemented? As far as I understand, they had only been declared before? They are still not implemented. The cmath header in libc++ provides wrappers around these functions for C++, but if they are not declared then the compiler throws an error. Now there is a linker error if they are used, but if they are not used then it works irrespective of whether you just include cmath or do undefined things like include math.h first. If this is right, are they usable on a recent CURRENT, built with gcc42 (system compiler), by ports which use gcc46 (not clang)? If not, are there any plans to implement these functions in the near future? I think they're one of the last bits of C99 / C11 that anyone actually cares about (C11 unicode support being another), so they'll probably get done at some point, but doing them correctly is nontrivial, except on platforms where double and long double are the same. The use of C99 extra long double functions would be of interest for example for programs like math/R, especially its upcoming releases. I would be very wary of anything that depends on long double. The C spec makes no constraints on the range of float, double, or long double, but by convention on most platforms float is an IEEE 754 32-bit floating point value and double is 64-bit. No such conventions apply to long doubles and I know of (widely deployed) platforms where they are 64 bits, 80 bits and 128 bits, respectively. I believe on PowerPC FreeBSD actually gets it wrong and uses 64 bits even though the platform ABI spec recommends using 128 bits. x86 uses 80-bit x87 values (although sizeof(long double) may be 16 because they have some padding bytes in memory). If your program is tolerant of a potential variation in precision between 64 bits and 128 bits, then it is probably tolerant of just using doubles. David___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 28.05.2012 10:41 (UTC+1), David Chisnall wrote: On 28 May 2012, at 05:35, Rainer Hurling wrote: First I should note that I am by no means an expert in C / C++ or C99 standards. So my comments are only on a level of someone who is using FreeBSD for scientific purposes like GIS and math applications and who is the maintainer of some scientific ports (like math/saga). Yesterday r236148 (Allow inclusion of libc++cmath to work after including math.h) was comitted to head, many thanks. Does this mean, that extra long double functions like acoshl, expm1l or log1pl are now really implemented? As far as I understand, they had only been declared before? They are still not implemented. Thecmath header in libc++ provides wrappers around these functions for C++, but if they are not declared then the compiler throws an error. Now there is a linker error if they are used, but if they are not used then it works irrespective of whether you just includecmath or do undefined things like includemath.h first. Ok, that's what I had supposed. Because the main difference between r236147 and r2136148 seems to be the define of _MATH_EXTRA_H_, the rest is more a type of binning? If this is right, are they usable on a recent CURRENT, built with gcc42 (system compiler), by ports which use gcc46 (not clang)? If not, are there any plans to implement these functions in the near future? I think they're one of the last bits of C99 / C11 that anyone actually cares about (C11 unicode support being another), so they'll probably get done at some point, but doing them correctly is nontrivial, except on platforms where double and long double are the same. Yes, I agree. These outstanding long double math functions (like log1pl) and better unicode support are really needed for some important third party projects like R or SAGA GIS. In the past I have read several times discussions about the correctness of long double functions on FreeBSD. Some drafts have been dismissed because of there inaccuracy in special cases. Also was discussed to get missing libm routines from NetBSD [1]. It appears as if we have to wait some more time ... The use of C99 extra long double functions would be of interest for example for programs like math/R, especially its upcoming releases. I would be very wary of anything that depends on long double. The C spec makes no constraints on the range of float, double, or long double, but by convention on most platforms float is an IEEE 754 32-bit floating point value and double is 64-bit. No such conventions apply to long doubles and I know of (widely deployed) platforms where they are 64 bits, 80 bits and 128 bits, respectively. I believe on PowerPC FreeBSD actually gets it wrong and uses 64 bits even though the platform ABI spec recommends using 128 bits. x86 uses 80-bit x87 values (although sizeof(long double) may be 16 because they have some padding bytes in memory). If your program is tolerant of a potential variation in precision between 64 bits and 128 bits, then it is probably tolerant of just using doubles. Yes, I think in most cases math/R is tolerant enough of just using doubles. But in the near future they plan to implement much more of the C99 stuff and their tolerance to offer workarounds for FreeBSD shrinks from release to release [2]. So these problems will increase :-( Many thanks for your fast und well explained answer, I really appreciate it. Rainer David [1] http://lists.freebsd.org/pipermail/freebsd-standards/2011-February/002119.html [2] https://stat.ethz.ch/pipermail/r-devel/2012-May/064128.html And some more places I found interesting about C99 and FreeBSD [3] http://www.freebsd.org/projects/c99/index.html [4] http://wiki.freebsd.org/MissingMathStuff [5] http://marc.info/?l=freebsd-standardsm=123158761305427 [6] http://lists.freebsd.org/pipermail/freebsd-hackers/2009-March/028030.html [7] http://www.koders.com/c/fid164FD2F50C9AAA63119A641875824455B3AE6B55.aspx?s=log1pl.c [8] http://www.mail-archive.com/bug-gnulib@gnu.org/msg26571.html [9] http://leaf.dragonflybsd.org/mailarchive/commits/2011-12/msg00190.html ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 28 May 2012, at 13:30, Rainer Hurling wrote: On 28.05.2012 10:41 (UTC+1), David Chisnall wrote: On 28 May 2012, at 05:35, Rainer Hurling wrote: Ok, that's what I had supposed. Because the main difference between r236147 and r2136148 seems to be the define of _MATH_EXTRA_H_, the rest is more a type of binning? Yes, it's just about making libc++'s cmath header compile, nothing more. Yes, I agree. These outstanding long double math functions (like log1pl) and better unicode support are really needed for some important third party projects like R or SAGA GIS. I very much doubt that anything is using the C11 unicode stuff yet, since no compiler or libc currently supports it... In the past I have read several times discussions about the correctness of long double functions on FreeBSD. Some drafts have been dismissed because of there inaccuracy in special cases. Also was discussed to get missing libm routines from NetBSD [1]. It appears as if we have to wait some more time ... I thought we did pull in some NetBSD libm stuff recently. Not sure what the status of that is, you'd need to check with das@. Yes, I think in most cases math/R is tolerant enough of just using doubles. But in the near future they plan to implement much more of the C99 stuff and their tolerance to offer workarounds for FreeBSD shrinks from release to release [2]. So these problems will increase :-( Reading that email, it seems that they would prefer a function that exists and returns the wrong result to one that does not exist. If this is really the attitude of the developers of R, then I shall make absolutely certain that I avoid using R for anything where I actually care about the results, and would strongly encourage everyone else to do the same. In general, (sane) people use the long double versions because they need the extra precision and care about the result. We could easily implement the long double versions now as toy versions that just called the double versions, but that would upset a lot of people who actually care about accuracy, who are the main target audience for these functions. David___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 28.05.2012 14:49 (UTC+1), David Chisnall wrote: On 28 May 2012, at 13:30, Rainer Hurling wrote: On 28.05.2012 10:41 (UTC+1), David Chisnall wrote: On 28 May 2012, at 05:35, Rainer Hurling wrote: Ok, that's what I had supposed. Because the main difference between r236147 and r2136148 seems to be the define of _MATH_EXTRA_H_, the rest is more a type of binning? Yes, it's just about making libc++'s cmath header compile, nothing more. I see, thanks. Yes, I agree. These outstanding long double math functions (like log1pl) and better unicode support are really needed for some important third party projects like R or SAGA GIS. I very much doubt that anything is using the C11 unicode stuff yet, since no compiler or libc currently supports it... Of course you are right with C11 unicode stuff. I thougt more about my actual problems to get the right charset conversions between different apps, i.e. qgis - wxgtk29 - saga gis. Or using Gnome apps (utf8) on windowmaker using ISO8859-15. But this is OT here, sorry. In the past I have read several times discussions about the correctness of long double functions on FreeBSD. Some drafts have been dismissed because of there inaccuracy in special cases. Also was discussed to get missing libm routines from NetBSD [1]. It appears as if we have to wait some more time ... I thought we did pull in some NetBSD libm stuff recently. Not sure what the status of that is, you'd need to check with das@. I am not aware of it and will have a look. But this did not implement the missing long double functions. Yes, I think in most cases math/R is tolerant enough of just using doubles. But in the near future they plan to implement much more of the C99 stuff and their tolerance to offer workarounds for FreeBSD shrinks from release to release [2]. So these problems will increase :-( Reading that email, it seems that they would prefer a function that exists and returns the wrong result to one that does not exist. If this is really the attitude of the developers of R, then I shall make absolutely certain that I avoid using R for anything where I actually care about the results, and would strongly encourage everyone else to do the same. This was a statement of just one (though not unimportant) person from the R Core team. I don't think that this is the only view of R Core developers. On the other hand he is the person, who actually did most of the stuff within R for years now to circumvent the problems running on FreeBSD. In general, (sane) people use the long double versions because they need the extra precision and care about the result. We could easily implement the long double versions now as toy versions that just called the double versions, but that would upset a lot of people who actually care about accuracy, who are the main target audience for these functions. It seems to be a general trend (outside of FreeBSD) to implement more and more stuff at the cost of quality. I am certain that for many FreeBSD users accuracy is more important than completeness, especially for scientists. Nevertheless this policy brings some problems in the real world ;-) Thanks again for your thoughts, Rainer David ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
One thing that could be done is to have a math/cephes port that adds the extra C99 math functions. This is already done in the math/sage port, using a rather clever patch due to Peter Jeremy, that applies to the cephes code. What it would do is to create a /usr/local/lib/libm.so that would provide the extra functions not currently included in /lib/libm.so, and then link in /lib/libm.so as well. It would also create its own /usr/local/include/math.h and /usr/local/include/complex.h as well. What do you guys think? Do you want someone to start experimenting with this idea? I could do it, but probably not for a little while. Stephen ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Mon, May 28, 2012 at 11:01:24AM -0500, Stephen Montgomery-Smith wrote: One thing that could be done is to have a math/cephes port that adds the extra C99 math functions. This is already done in the math/sage port, using a rather clever patch due to Peter Jeremy, that applies to the cephes code. What it would do is to create a /usr/local/lib/libm.so that would provide the extra functions not currently included in /lib/libm.so, and then link in /lib/libm.so as well. It would also create its own /usr/local/include/math.h and /usr/local/include/complex.h as well. What do you guys think? Do you want someone to start experimenting with this idea? I could do it, but probably not for a little while. This is a horrible, horrible, horrible idea. Have you looked at the cephes code, particularly the complex.h functions? -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 05/28/2012 03:31 PM, Steve Kargl wrote: On Mon, May 28, 2012 at 11:01:24AM -0500, Stephen Montgomery-Smith wrote: One thing that could be done is to have a math/cephes port that adds the extra C99 math functions. This is already done in the math/sage port, using a rather clever patch due to Peter Jeremy, that applies to the cephes code. What it would do is to create a /usr/local/lib/libm.so that would provide the extra functions not currently included in /lib/libm.so, and then link in /lib/libm.so as well. It would also create its own /usr/local/include/math.h and /usr/local/include/complex.h as well. What do you guys think? Do you want someone to start experimenting with this idea? I could do it, but probably not for a little while. This is a horrible, horrible, horrible idea. Have you looked at the cephes code, particularly the complex.h functions? I have only taken a very cursory look. What should I specifically look for in seeing that the code is bad? ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 2012-May-28 11:01:24 -0500, Stephen Montgomery-Smith step...@missouri.edu wrote: One thing that could be done is to have a math/cephes port that adds the extra C99 math functions. This is already done in the math/sage port, using a rather clever patch due to Peter Jeremy, that applies to the cephes code. What it would do is to create a /usr/local/lib/libm.so that would provide the extra functions not currently included in /lib/libm.so, and then link in /lib/libm.so as well. It would also create its own /usr/local/include/math.h and /usr/local/include/complex.h as well. Basically, as long as the compiler searches /usr/local/{include,lib} before the base include/lib then math.h, complex.h and -lm give the application a complete C99 math implementation by using base functions where they exist and cephes functions where they don't. The patch I wrote for sage can be found at http://trac.sagemath.org/sage_trac/ticket/9543 If there's any interest, I could produce a port for this. Another option would be to import cephes into base and use it to provide the missing C99 functions. Cephes includes copyright notices but the closest I can find to a license is: Some software in this archive may be from the book _Methods and Programs for Mathematical Functions_ (Prentice-Hall or Simon Schuster International, 1989) or from the Cephes Mathematical Library, a commercial product. In either event, it is copyrighted by the author. What you see here may be used freely but it comes with no support or guarantee. -- Peter Jeremy pgpYmCz2gMd3i.pgp Description: PGP signature
Re: Use of C99 extra long double math functions after r236148
On 2012-May-28 13:31:59 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: On Mon, May 28, 2012 at 11:01:24AM -0500, Stephen Montgomery-Smith wrote: One thing that could be done is to have a math/cephes port that adds the extra C99 math functions. This is already done in the math/sage port, using a rather clever patch due to Peter Jeremy, that applies to the cephes code. ... This is a horrible, horrible, horrible idea. Have you looked at the cephes code, particularly the complex.h functions? The cephes code is somewhat a mess layout-wise. Algorithmetically, it seems somewhat variable - some functions are implemented (hopefully correctly) using semi-numerical techniques, whereas others just use mathematical identities which will result in precision loss - though most of the functions include accuracy information. I agree it would be far preferable to have a properly validated C99 libm with all functions having maximum errors of a no more than a few LSB over their complete domain, as well as correct support for signed zeroes, infinities and signalling and non-signalling NaNs but that is a non-trivial undertaking. In the interim, how should FreeBSD handle apps that want a C99 libm? 1) Fail to build them 2) Provide possibly imperfect fallbacks for the unimplemented bits. If someone (I don't have the expertise) wants to identify the cephes functions that are sub-standard, we can include link-time warnings (as done for eg gets(3)) when they are used. -- Peter Jeremy pgpcG5SKNkFm9.pgp Description: PGP signature
Re: Use of C99 extra long double math functions after r236148
On Mon, May 28, 2012 at 04:19:22PM -0500, Stephen Montgomery-Smith wrote: On 05/28/2012 03:31 PM, Steve Kargl wrote: On Mon, May 28, 2012 at 11:01:24AM -0500, Stephen Montgomery-Smith wrote: One thing that could be done is to have a math/cephes port that adds the extra C99 math functions. This is already done in the math/sage port, using a rather clever patch due to Peter Jeremy, that applies to the cephes code. What it would do is to create a /usr/local/lib/libm.so that would provide the extra functions not currently included in /lib/libm.so, and then link in /lib/libm.so as well. It would also create its own /usr/local/include/math.h and /usr/local/include/complex.h as well. What do you guys think? Do you want someone to start experimenting with this idea? I could do it, but probably not for a little while. This is a horrible, horrible, horrible idea. Have you looked at the cephes code, particularly the complex.h functions? I have only taken a very cursory look. What should I specifically look for in seeing that the code is bad? Well, to start with, the extra C99 math functions that are missing in libm include a few for the long double type and many (if not most) of the complex functions for any type. Cephes at best will give you float, double, and ld80, but not ld128 versions of the functions. Then, there is fun little fact that neither (base) gcc nor clang nor gcc less than 4.6 does complex arithematic correctly; so, one needs to jump through hoops to get the correct answer (See Annex G in n1256.pdf). One item of particular importance in Annex G is the behavior of the functions for Re(z) and/or Im(z) being +-0, +-Inf, and NaN. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Tue, May 29, 2012 at 08:04:36AM +1000, Peter Jeremy wrote: On 2012-May-28 13:31:59 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: On Mon, May 28, 2012 at 11:01:24AM -0500, Stephen Montgomery-Smith wrote: One thing that could be done is to have a math/cephes port that adds the extra C99 math functions. This is already done in the math/sage port, using a rather clever patch due to Peter Jeremy, that applies to the cephes code. ... This is a horrible, horrible, horrible idea. Have you looked at the cephes code, particularly the complex.h functions? The cephes code is somewhat a mess layout-wise. Algorithmetically, it seems somewhat variable - some functions are implemented (hopefully correctly) using semi-numerical techniques, whereas others just use mathematical identities which will result in precision loss - though most of the functions include accuracy information. I agree it would be far preferable to have a properly validated C99 libm with all functions having maximum errors of a no more than a few LSB over their complete domain, as well as correct support for signed zeroes, infinities and signalling and non-signalling NaNs but that is a non-trivial undertaking. In the interim, how should FreeBSD handle apps that want a C99 libm? 1) Fail to build them 2) Provide possibly imperfect fallbacks for the unimplemented bits. If someone (I don't have the expertise) wants to identify the cephes functions that are sub-standard, we can include link-time warnings (as done for eg gets(3)) when they are used. Given that cephes was written years before C99 was even conceived, I suspect all functions are sub-standard. For example, AFAIK, none of the long double functions are appropriate for any platform that has an 128-bit long double; as cephes was written for an Intel 80-bit format. If portmgr or a port maintainer wants to use a library with untested implementations of missing libm functions, please do not put it into /usr/local/lib and call it libm. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Tue, May 29, 2012 at 07:05:07AM +1000, Peter Jeremy wrote: On 2012-May-28 11:01:24 -0500, Stephen Montgomery-Smith step...@missouri.edu wrote: One thing that could be done is to have a math/cephes port that adds the extra C99 math functions. This is already done in the math/sage port, using a rather clever patch due to Peter Jeremy, that applies to the cephes code. What it would do is to create a /usr/local/lib/libm.so that would provide the extra functions not currently included in /lib/libm.so, and then link in /lib/libm.so as well. It would also create its own /usr/local/include/math.h and /usr/local/include/complex.h as well. Basically, as long as the compiler searches /usr/local/{include,lib} before the base include/lib then math.h, complex.h and -lm give the application a complete C99 math implementation by using base functions where they exist and cephes functions where they don't. The patch I wrote for sage can be found at http://trac.sagemath.org/sage_trac/ticket/9543 If there's any interest, I could produce a port for this. Another option would be to import cephes into base and use it to provide the missing C99 functions. Cephes includes copyright notices but the closest I can find to a license is: Some software in this archive may be from the book _Methods and Programs for Mathematical Functions_ (Prentice-Hall or Simon Schuster International, 1989) or from the Cephes Mathematical Library, a commercial product. In either event, it is copyrighted by the author. What you see here may be used freely but it comes with no support or guarantee. Please talk to das@ (although I believe he's finishing up his dissertation). I recall that he's stated that he looked into using cephes, and concluded that it is not suitable for libm. Note there is also http://www.freebsd.org/cgi/query-pr.cgi?pr=kern/147599 which I've also objected to importing into libm. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 05/28/2012 05:17 PM, Steve Kargl wrote: On Mon, May 28, 2012 at 04:19:22PM -0500, Stephen Montgomery-Smith wrote: On 05/28/2012 03:31 PM, Steve Kargl wrote: On Mon, May 28, 2012 at 11:01:24AM -0500, Stephen Montgomery-Smith wrote: One thing that could be done is to have a math/cephes port that adds the extra C99 math functions. This is already done in the math/sage port, using a rather clever patch due to Peter Jeremy, that applies to the cephes code. What it would do is to create a /usr/local/lib/libm.so that would provide the extra functions not currently included in /lib/libm.so, and then link in /lib/libm.so as well. It would also create its own /usr/local/include/math.h and /usr/local/include/complex.h as well. What do you guys think? Do you want someone to start experimenting with this idea? I could do it, but probably not for a little while. This is a horrible, horrible, horrible idea. Have you looked at the cephes code, particularly the complex.h functions? I have only taken a very cursory look. What should I specifically look for in seeing that the code is bad? Well, to start with, the extra C99 math functions that are missing in libm include a few for the long double type and many (if not most) of the complex functions for any type. Cephes at best will give you float, double, and ld80, but not ld128 versions of the functions. Then, there is fun little fact that neither (base) gcc nor clang nor gcc less than 4.6 does complex arithematic correctly; so, one needs to jump through hoops to get the correct answer (See Annex G in n1256.pdf). One item of particular importance in Annex G is the behavior of the functions for Re(z) and/or Im(z) being +-0, +-Inf, and NaN. IMHO these problems are definitely not bad enough to avoid making a math/cephes port. I for one would definitely like to have some kind of implementation of ccosh, even if it gets a few things wrong when it is fed Nan or I*Inf or such like as its input. I mean, if clang or gcc46 doesn't even get this right, how can we expect better from libm? Also, a really nice thing about Jeremy's patch is that it automatically detects which functions are already included in the base libm, and only adds functions not already in libm. And ld80 is better than nothing if ld128 isn't available. I would avoid cephes only if it got some of the answers horribly wrong (as in erf(100) being something like -14232 as the openoffice implementation of the erf function used to report). I say, lets go ahead and add a math/cephes port. By the way, do you have an opinion on the complex functions used in Linux? (I tried reading the glibc code, but I could only find chains of macros and functions calling each other, and I could never find where the actual work was performed.) ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Mon, May 28, 2012 at 06:03:37PM -0500, Stephen Montgomery-Smith wrote: On 05/28/2012 05:17 PM, Steve Kargl wrote: On Mon, May 28, 2012 at 04:19:22PM -0500, Stephen Montgomery-Smith wrote: On 05/28/2012 03:31 PM, Steve Kargl wrote: On Mon, May 28, 2012 at 11:01:24AM -0500, Stephen Montgomery-Smith wrote: One thing that could be done is to have a math/cephes port that adds the extra C99 math functions. This is already done in the math/sage port, using a rather clever patch due to Peter Jeremy, that applies to the cephes code. What it would do is to create a /usr/local/lib/libm.so that would provide the extra functions not currently included in /lib/libm.so, and then link in /lib/libm.so as well. It would also create its own /usr/local/include/math.h and /usr/local/include/complex.h as well. What do you guys think? Do you want someone to start experimenting with this idea? I could do it, but probably not for a little while. This is a horrible, horrible, horrible idea. Have you looked at the cephes code, particularly the complex.h functions? I have only taken a very cursory look. What should I specifically look for in seeing that the code is bad? Well, to start with, the extra C99 math functions that are missing in libm include a few for the long double type and many (if not most) of the complex functions for any type. Cephes at best will give you float, double, and ld80, but not ld128 versions of the functions. Then, there is fun little fact that neither (base) gcc nor clang nor gcc less than 4.6 does complex arithematic correctly; so, one needs to jump through hoops to get the correct answer (See Annex G in n1256.pdf). One item of particular importance in Annex G is the behavior of the functions for Re(z) and/or Im(z) being +-0, +-Inf, and NaN. IMHO these problems are definitely not bad enough to avoid making a math/cephes port. I for one would definitely like to have some kind of implementation of ccosh, even if it gets a few things wrong when it is fed Nan or I*Inf or such like as its input. I mean, if clang or gcc46 doesn't even get this right, how can we expect better from libm? Search the llvm and gcc bug databases, I'm the person that informed both that they can't deal with complex arithematic correctly. Interesting that you mention ccosh. From clog.c in http://www.netlib.org/cephes/cmath.tgz void ccosh (z, w) cmplx *z, *w; { double x, y; x = z-r; y = z-i; w-r = cosh (x) * cos (y); w-i = sinh (x) * sin (y); } From clog.c in http://www.netlib.org/cephes/c9x-complex double complex ccosh (z) double complex z; { double complex w; double x, y; x = creal(z); y = cimag(z); w = cosh (x) * cos (y) + (sinh (x) * sin (y)) * I; return (w); } See math_private.h about the above. And, finally, : /*- * Copyright (c) 2005 Bruce D. Evans and Steven G. Kargl * All rights reserved. * * Redistribution and use in source and binary forms, with or without * modification, are permitted provided that the following conditions * are met: * 1. Redistributions of source code must retain the above copyright *notice unmodified, this list of conditions, and the following *disclaimer. * 2. Redistributions in binary form must reproduce the above copyright *notice, this list of conditions and the following disclaimer in the *documentation and/or other materials provided with the distribution. * * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT, * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. */ /* * Hyperbolic cosine of a complex argument z = x + i y. * * cosh(z) = cosh(x+iy) * = cosh(x) cos(y) + i sinh(x) sin(y). * * Exceptional values are noted in the comments within the source code. * These values and the return value were taken from n1124.pdf. */ #include sys/cdefs.h __FBSDID($FreeBSD: head/lib/msun/src/s_ccosh.c 226599 2011-10-21 06:29:32Z das $); #include complex.h #include math.h #include math_private.h static const double huge = 0x1p1023; double complex ccosh(double complex z) { double x, y, h; int32_t hx, hy, ix, iy, lx, ly; x = creal(z); y = cimag(z); EXTRACT_WORDS(hx, lx, x); EXTRACT_WORDS(hy, ly, y); ix = 0x7fff hx; iy = 0x7fff hy; /* Handle the
Re: Use of C99 extra long double math functions after r236148
On 05/28/2012 06:30 PM, Steve Kargl wrote: From clog.c in http://www.netlib.org/cephes/c9x-complex double complex ccosh (z) double complex z; { double complex w; double x, y; x = creal(z); y = cimag(z); w = cosh (x) * cos (y) + (sinh (x) * sin (y)) * I; return (w); } See math_private.h about the above. I looked in math_private.h - I presume you meant lib/msun/src/math_private.h. I wasn't able to find anything about ccosh there. I think that for a rough and ready ccosh, this is high enough quality for a math/cephes port. I do agree that it might not be high enough quality to make FreeBSD base. (Although I don't think I have ever been in a situation where I would have been tripped up by a transcendental function that responded incorrectly to exceptional input.) And, finally, Yes, it is very nice. Who's writing the code to test the implementations? That is better much the problem. Without testing, one might get an implementation that appears to work until it doesn't! It took me 3+ years to get sqrtl() into libm, but bde and das (and myself) wanted to make sure the code worked. Fair enough if we are talking about the base system. I haven't looked at glibc code in years, because I hack on libm when I can. I do not want to run into questions about whether my code is tainted by the gpl. They had similar lists of exceptions. ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On Mon, May 28, 2012 at 06:44:42PM -0500, Stephen Montgomery-Smith wrote: On 05/28/2012 06:30 PM, Steve Kargl wrote: From clog.c in http://www.netlib.org/cephes/c9x-complex double complex ccosh (z) double complex z; { double complex w; double x, y; x = creal(z); y = cimag(z); w = cosh (x) * cos (y) + (sinh (x) * sin (y)) * I; return (w); } See math_private.h about the above. I looked in math_private.h - I presume you meant lib/msun/src/math_private.h. I wasn't able to find anything about ccosh there. I think that for a rough and ready ccosh, this is high enough quality for a math/cephes port. That's the problem. It is not acceptable (imo). The comment in math_private.h that is relevant is /* * Inline functions that can be used to construct complex values. * * The C99 standard intends x+I*y to be used for this, but x+I*y is * currently unusable in general since gcc introduces many overflow, * underflow, sign and efficiency bugs by rewriting I*y as * (0.0+I)*(y+0.0*I) and laboriously computing the full complex product. * In particular, I*Inf is corrupted to NaN+I*Inf, and I*-0 is corrupted * to -0.0+I*0.0. */ Those wrong +-0 mean you may up end up on the worng riemann sheet, and that NaN propagates. -- Steve ___ freebsd-current@freebsd.org mailing list http://lists.freebsd.org/mailman/listinfo/freebsd-current To unsubscribe, send any mail to freebsd-current-unsubscr...@freebsd.org
Re: Use of C99 extra long double math functions after r236148
On 2012-May-28 15:54:06 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: Given that cephes was written years before C99 was even conceived, I suspect all functions are sub-standard. Well, most of cephes was written before C99. The C99 parts of cephes were written to turn it into a complete C99 implementation. For example, AFAIK, none of the long double functions are appropriate for any platform that has an 128-bit long double; as cephes was written for an Intel 80-bit format. FreeBSD currently supports: 64-bit long doubles on ARM, MIPS and PowerPC; 80-bit long doubles on amd64, i386 and iA64; 128-bit long doubles on SPARC. The lack of LD128 in cephes therefore only affects one (not widely used) platform. The lack of even de facto standards for long double mean that any applications wanting to use them already need to cope with at least a 2:1 precision range. If portmgr or a port maintainer wants to use a library with untested implementations of missing libm functions, please do not put it into /usr/local/lib and call it libm. There some test code in cephes. Can you point me to a suitable test suite for LD80 and LD128? The reason for calling it libm is to avoid having to hack every consumer to add an additional library. On 2012-May-28 16:30:35 -0700, Steve Kargl s...@troutmask.apl.washington.edu wrote: Who's writing the code to test the implementations? That is better much the problem. Without testing, one might get an implementation that appears to work until it doesn't! That is equally true of the rest of FreeBSD. The list of open PRs suggests that FreeBSD still has a fair way to go before reaching perfection. And, most of this thread has been about using this code in ports - where the bar is much lower. Who is writing the code to test all the other ports? What is so special about this particular proposed port that it needs to come with solid-gold credentials? It took me 3+ years to get sqrtl() into libm, but bde and das (and myself) wanted to make sure the code worked. Last time I checked (a couple of years ago), FreeBSD was missing 65 C99 libm functions. At 3 years per function, we should have C99 support available early in the 23rd century - which may be a bit late. On 2012-May-28 22:03:43 -0500, Stephen Montgomery-Smith step...@missouri.edu wrote: 1. By being so picky about being so precise, FreeBSD is behind the time line in rolling out a usable set of C99 functions. And at the current rate, we'll all be long dead before they are available. Whilst I'd far prefer to have a properly verifed library function, I think we are better off with an implementation that has some caveats regarding edge-case behaviour than having nothing. In the end, I do think it is good to ultimately settle on good C99 compliant code. But having something intermediate that mostly works is better than nothing. Especially if it exists only in the ports, and not in the base code. I agree with this sentiment. What do people do on other free OSs? Does a tested open source C99 libm exist anywhere? glibc implements cpow(x,y) as cexp(y*clog(x)) and cephes does better than that. Is FreeBSD wasting its time writing correct C99 code because all the libm consumers expect no better than what glibc offers? I agree that writing correct libm functions is hard. I think a lot of the problem is that it's a mix of lots of boilerplate code testing for special conditions and edge cases that is boring to write and fiddly to get right, together with a kernel that is a pile of polynomial evaluations full of magic numbers that needs specialist skills to write. If we could get someone with the relevant skills to formally list all the special conditions edge cases for each function, it should be possible to generate both the library C code and test cases from that - which would remove a lot of the tedium. -- Peter Jeremy pgpUnZGDcc79l.pgp Description: PGP signature