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.orig    2012-07-25 09:38:31.000000000 -0700
> +++ a.c 2012-07-25 09:40:36.000000000 -0700
> _at__at_ -1,5 +1,6 _at__at_
>   #include <stdio.h>
>   #include <math.h>
> +#include <float.h>
>
>   int main(void)
>   {
> _at__at_ -9,8 +10,8 _at__at_
>     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.