Problems to calculate sin

[Kenneth P. Turvey]

On Wed, 05 Mar 2008 01:09:56 +0000, Roedy Green wrote:

[Snip]

the Intel FP instruction set has a sine-computing instruction. It
works inside with polynomial approximations. Any error is Intel's
doing.

[Snip]

I read somewhere that Java doesn't use this sine-computing instruction
since it doesn't meet the accuracy requirements guaranteed by the class.
This was quoted as the reason that Java performs slower on transcendental
functions than other languages on the platform.

Is this information out of date?

 

[Thomas.a.mcglynn]

(e-mail address removed) wrote:

...


How many bytes per interval?

Patricia

For a cubic fit presumably one needs 4 values per interval or 32 bytes
for double coefficients. The expansion isn't around 0 so the even and
odd terms both show. So about 300kB altogether. Too big to want to
use it without a good reason, but if, contrary to our belief, the
computation of the sine were faster using the table, easily
accommodated within typical programs, and a bit smaller than the
typical cache size. Presumably one would need to compute a lot of
sines to accommodate the setup overhead though.

Tom

 

[Patricia Shanahan]

For a cubic fit presumably one needs 4 values per interval or 32 bytes
for double coefficients. The expansion isn't around 0 so the even and
odd terms both show. So about 300kB altogether. Too big to want to
use it without a good reason, but if, contrary to our belief, the
computation of the sine were faster using the table, easily
accommodated within typical programs, and a bit smaller than the
typical cache size. Presumably one would need to compute a lot of
sines to accommodate the setup overhead though.

Tom

I think at the moment the polynomial approximation approach wins,
because the processor has the greatest freedom to reorder and schedule
arithmetic involving only a few constants and one parameter, with no
memory loads.

As caches get bigger, the economics may shift again.

It's the sort of thing where if I were doing it professionally I might
maintain implementations both ways, and compare them on each processor
generation to see shifts in the relative performance.

Patricia

 

[Thomas.a.mcglynn]

....

I think at the moment the polynomial approximation approach wins,
because the processor has the greatest freedom to reorder and schedule
arithmetic involving only a few constants and one parameter, with no
memory loads.

As caches get bigger, the economics may shift again.

It's the sort of thing where if I were doing it professionally I might
maintain implementations both ways, and compare them on each processor
generation to see shifts in the relative performance.

Patricia

I'm sure your right about the current situation. However I suspect I
overestimated the needed size for a cache somewhat. I was naively
thinking of 10,000 independent cubic fits, but it's more reasonable
and efficient to use the values from several points near where the
interpolation is to be done. In that case we only need a single
coefficient for each point (presumably the function value) so we can
get by with a mere ~100kB of table. I don't think this affects your
conclusion though.

Regards,
Tom

 

[Mark Thornton]

On Wed, 05 Mar 2008 01:09:56 +0000, Roedy Green wrote:

[Snip]

the Intel FP instruction set has a sine-computing instruction. It
works inside with polynomial approximations. Any error is Intel's
doing.

[Snip]

I read somewhere that Java doesn't use this sine-computing instruction
since it doesn't meet the accuracy requirements guaranteed by the class.
This was quoted as the reason that Java performs slower on transcendental
functions than other languages on the platform.

It does use the FSIN function, but only after reducing the argument to
+-PI/4 (I believe). The accuracy problem with the FSIN function is the
way Intel do argument reduction using a 66bit approximation of PI. The
Java specification requires many more bits of PI in some cases.

Mark Thornton

 

[Mark Thornton]

I'm sure your right about the current situation. However I suspect I
overestimated the needed size for a cache somewhat. I was naively
thinking of 10,000 independent cubic fits, but it's more reasonable
and efficient to use the values from several points near where the
interpolation is to be done. In that case we only need a single
coefficient for each point (presumably the function value) so we can
get by with a mere ~100kB of table. I don't think this affects your
conclusion though.

Regards,
Tom

It may be worth noting that the Level one data caches are often much
smaller than this (32KB/core on my machine).

Mark Thornton

 

[Thomas.a.mcglynn]

I suppose in some future chip it will have special in-parallel checks
for 45 degrees, 90 degrees, 0 degrees, 30 degrees to get as perfect as
possible results.

Hmmm.... Do you mean that you would like something like

Math.sin(Math.PI)

to return 0? That seems very dangerous.

I could see having special checking on new degree based functions
though: E.g., [untested]

double sind(double x) {
double specialValue[] {
0,0.5,Math.sqrt(3)/2, 1, ...};

x = x % 360;
if (x < 0) {
x += 360;
}

if ( (x mod 30) == 0} {
return specialValue[x/30];
} else {
return Math.sin(x*Math.PI/180);
}
}

Regards,
Tom McGlynn

 

[Roedy Green]

I read somewhere that Java doesn't use this sine-computing instruction
since it doesn't meet the accuracy requirements guaranteed by the class.
This was quoted as the reason that Java performs slower on transcendental
functions than other languages on the platform.

Is this information out of date?

I don't know. The code would in a native class. The chip will not have
changed, so likely that information would be stable.

 

[Andreas Leitgeb]

On 05 Mar 2008 09:41:15 GMT, "Kenneth P. Turvey"

I don't know. The code would in a native class. The chip will not have
changed, so likely that information would be stable.

My guess was, that if one needed math more exact than native
processor's arithmetics, he would use "strictfp", so while I
don't know for sure, I wouldn't bet on abovementioned information.

 

[Patricia Shanahan]

My guess was, that if one needed math more exact than native
processor's arithmetics, he would use "strictfp", so while I
don't know for sure, I wouldn't bet on abovementioned information.

The minimum precision requirements in the Math.sin API documentation are
not dependent on strictfp.

Patricia

 

[Tom McGlynn]

Andreas Leitgeb wrote: ....


The minimum precision requirements in the Math.sin API documentation are
not dependent on strictfp.

Perhaps Andreas may be thinking of the StrictMath class rather than
strictfp? This class duplicates the functionality of Math but
requires that the results be identical to the fdlibm results. [I
think somewhere else in this thread I also confused this with
strictfp].

It's a bit counterintuitive, but I believe the general effect of
either a strictfp block or use of the StrictMath class is a less
accurate calculation. They are all about consistency, not accuracy.

Regards,
Tom McGlynn

 

[Patricia Shanahan]

Andreas Leitgeb wrote: ...

The minimum precision requirements in the Math.sin API documentation are
not dependent on strictfp.

Perhaps Andreas may be thinking of the StrictMath class rather than
strictfp? This class duplicates the functionality of Math but
requires that the results be identical to the fdlibm results. [I
think somewhere else in this thread I also confused this with
strictfp].

It's a bit counterintuitive, but I believe the general effect of
either a strictfp block or use of the StrictMath class is a less
accurate calculation. They are all about consistency, not accuracy.

I think StrictMath can go either way. There are generally two
representable numbers that could be a valid Math.sin result for a given
angle, because each is within one ULP of the infinite precision answer.
Math.sin can return either. StrictMath.sin must return the one fdlibm
would pick, which may or may not be the closer of the two.

Patricia

 

[Steve70]

And why Turbo Pascal calculate sin(30*PI/180)=0.5?

 

[Lew]

And why Turbo Pascal calculate sin(30*PI/180)=0.5?

Does it? Does it really?

You can make the Java result of Math.sin( 30.0 * Math.PI / 180.0 ) look like
it returns 0.5, too. The significant thing is that the result of that
calculation has guaranteed results, precision and accuracy in Java. Your
Pascal or C implementations will vary in their internal precision, accuracy
and output.

As pointed out quite frequently in these newsgroups, floating point
expressions in a binary device will generally be inaccurate. The science of
numerical analysis allows limits on that inaccuracy, and yields algorithms to
minimize it. Programmers need to have some knowledge of that science in order
to understand, and more importantly, control what their code does.

Note the definition of Java's Math.PI:

The double value that is closer than any other to pi,
the ratio of the circumference of a circle to its diameter.

If PI cannot be exactly represented, then it is not reasonable to expect that
the result of a calculation on that approximation involving division by 180.0,
whose result generally will also be inaccurate, much less the result of an
approximation of a transcendental function of that calculation, will be
exactly representable in the binary floating-point format, much less exactly
accurate.

 

[Mark Thornton]

And why Turbo Pascal calculate sin(30*PI/180)=0.5?

And what result does it print for sin(30*PI/180)-0.5

Did you verify that the result was exactly 0.5 or merely converted to
string as 0.5.

 

[Andrew Thompson]

Math.sin(30*Math.PI/180)=0.49999999999999994

What's the problem?

Digital representation of numbers. Please
research IEEE 754.

I expect turbo Pascal was doing internal rounding
(knowing Borland* - suspecting at time of output).

* I cannot imagine that Turbo Pascal was *not*
using IEEE 754 internally - I wrote a few trivial
number crunching progs. using turbo Pascal and
found the results quite plausible (they confirmed
basic theories of physics I was trying to disprove).

 

[Andrew Thompson]

...I wrote a few trivial
number crunching progs. using turbo Pascal and
found the results quite plausible (they confirmed ..

..to 10+ decimal places..

 

[Andreas Leitgeb]

And what result does it print for sin(30*PI/180)-0.5
Did you verify that the result was exactly 0.5 or merely converted to
string as 0.5.

I do find it believeable that TurboPascal produces a more exact
result, since afaik it not only uses the processors "double extended"
floating point arithmetics(*), but also has 10bytes-datatype "extended"
to store such values rather than only the 8byte "double".

So, I'd expect even "sin(30*PI/180)-0.5" to be significantly
smaller (that is: nearer to zero) than in java, with or without
strictfp/StrictMath.

(*): From Wikipedia: IEEE 754 specifies four formats for representing
floating-point values: single-precision (32-bit), double-precision
(64-bit), single-extended precision (≥ 43-bit, not commonly used) and
double-extended precision (≥ 79-bit, usually implemented with 80 bits).


PS: don't have any TurboPascal here at hand. Only know it from back
in the MsDos-era.

 

[John W. Kennedy]

As pointed out quite frequently in these newsgroups, floating point
expressions in a binary device will generally be inaccurate.

In this case, binary has nothing to do with it. With binary /or/
decimal, no computer offers infinite precision. In fact, all other
things being equal, binary is better; decimal is only to be preferred in
a computational context with decimal quantization (i.e., finance).

 

[Lew]

In this case, binary has nothing to do with it. With binary /or/
decimal, no computer offers infinite precision. In fact, all other
things being equal, binary is better; decimal is only to be preferred in
a computational context with decimal quantization (i.e., finance).

Oo-kaay.

Binariness is not relevant in the abstract, but since this community uses Java
on binary devices it is relevant in practice. All that stuff about decimal,
not that it'd matter if you'd discussed duodecimal or sexagintadecimal, is by
the wayside and doesn't alter the salient point a jot.

The point is that fixed-precision computations will be inaccurate for real
numbers, and error management therein is a vital concern.