Date: Wed, 25 Jul 2012 19:33:07 +0200 From: Rainer Hurling <rhurlin@gwdg.de> To: Steve Kargl <sgk@troutmask.apl.washington.edu> Cc: David Schultz <das@FreeBSD.org>, freebsd-current@FreeBSD.org, Bruce Evans <bde@FreeBSD.org> Subject: Re: Use of C99 extra long double math functions after r236148 Message-ID: <50102DD3.6010700@gwdg.de> In-Reply-To: <20120725170004.GA72338@troutmask.apl.washington.edu> References: <20120708124047.GA44061@zim.MIT.EDU> <210816F0-7ED7-4481-ABFF-C94A700A3EA0@bsdimp.com> <20120708233624.GA53462@troutmask.apl.washington.edu> <4FFBF16D.2030007@gwdg.de> <2A1DE516-ABB4-49D7-8C3D-2C4DA2D9FCF5@bsdimp.com> <4FFC412B.4090202@gwdg.de> <20120710151115.GA56950@zim.MIT.EDU> <4FFC5E5D.8000502@gwdg.de> <20120710225801.GB58778@zim.MIT.EDU> <50101EDE.6030509@gwdg.de> <20120725170004.GA72338@troutmask.apl.washington.edu>
index | next in thread | previous in thread | raw e-mail
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.orig 2012-07-25 09:38:31.000000000 -0700
> +++ a.c 2012-07-25 09:40:36.000000000 -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.000000) is
>> 7.38905609893065040
>> expl(2.000000) 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.000000000000000000e+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.
help
Want to link to this message? Use this URL: <https://mail-archive.FreeBSD.org/cgi/mid.cgi?50102DD3.6010700>
