Bessel function yn appears only single-precision on IA64

[BSSG]

When I call the following Bessel function in a 'c' program on a PA-RISC and then a IA-64 box:

yn( 2, 0.98726577450859631657 )

I get the following:

PA-32 => -1.68339139976380280000
IA-64 => -1.68339139819686622879

which differ starting at the 9th decimal. So is the IA-64 call only running float precision and not double precision?

Thanks.

--
Bob

[Gregory Fruth]

The "0000" at the end of the PA-RISC
seems to indicate that those digits
are beyond PA's precision. IA-64
supports the Intel 80-bit floating type;
perhaps libm uses it internally. For
comparison, running your program on
IA-32 (Win2k/Cygwin/gcc) gives
-1.68339139976380280928, which has the same
# of digits as your IA-64 result but which
is closer in value to your PA result.

Also, perhaps PHSS_32066 (Math Library
Cumulative Patch) will help.

[BSSG]

Okay but I believe they should match out to about the 15th or 16th decimal for double-precision, even if the last digits are trucated or whatever. Based on the response to a ticket I went ahead and opened this morning, I did try applying that patch and I also just tried compiling with the -fpeval=double option, but had the same outcome. Anyway I'll see what they come up with. Thanks.

[A. Clay Stephenson]

One thing I notice immediately is that you are passing in arguments with far more than double-precision granularity and are expecting meaningful results. The yn() function itself may be fine but the problem could be in the string to double conversion in the compiler or in the printf output conversion. In any event, your argument for yn() has too many significant figures unless you are working in quad precision.

[BSSG]

HP's C manual indicates the input will be converted to the nearest represented value appropriate for the data range. I've done the same type of thing with most of the other math routines between the IA-64 and the PA_RISC and never had a problem. But if that is the case, I'd expect the same conversion on both platforms, so wouldn't the result be the same even if there is a problem with the string to number translation?

--
Bob

[Stephen Keane]

Try a simple test like...

#include

int main()
{
double d = 1.2345678910111213141516;
float f = 1.2345678910111213141516;

printf("%.20f\n", d);
printf("%.20f\n", f);

return(0);
}

When I run compile (using gcc) and run (PA-RISC 11.11) this, I get

1.23456789101112130000
1.23456788063049320000

Which indicates I'm getting around 16 digits of precision for double and around 8 for float.

[BSSG]

Yup, so that's why I'm wondering if the IA-64 version of the yn Bessel function is only returning single precision, even though the man page shows it as double precision. (There's also a single precision version listed as ynf.)

BTW, I also get similar problems with the y1 function.

--
Bob

[Gregory Fruth]

Running the attached program on PA 2.0 and
IA-64 shows that they both get the same
internal bit patterns for 0.987265... (so
it's not a compiler parsing or type
conversion problem), and the value
generated by yn() is indeed a different
bit pattern (so it's not a printf()
problem). Curious.

[Dennis Handly]

>which differ starting at the 9th decimal. So is the IA-64 call only running float precision and not double precision?

If it was float, it would differ in the 7th decimal place.

It looks like you found a known problem:

PHSS_33276: s700_800 11.23 Math Library Cumulative Patch

JAGaf55860: Bessel functions have precision problems

The corrected values:

-1.68339139976380275131 long double

-1.68339139976380280928 double
-1.68339145183563232422 float

>IA-64 supports the Intel 80-bit floating type; perhaps libm uses it internally.

Naturally, the Integrity Math Lib was rewritten to be faster and more accurate.

>perhaps PHSS_32066 (Math Library Cumulative Patch) will help.

Unfortunately not until PHSS_33276.

Most math function are accurate to .5 ULP.  But Bessel and a few other functions don't meet this:

http://www.hp.com/go/fp#5.5