Quote:
Originally Posted by Thomas Richter
For 3.1.4, this did work, but this is due to the "round to zero" policy. That is, "round to zero in high precision" followed by "round to zero in low precision" is identical to "round to zero in low precision".
|
Ok, so we have either rounding issues by the round-to-zero policy, or we have rounding issues when using mathieeesingbas and mathieeedoubbas in parallel. Sigh.
Quote:
|
Unfortunately, "round to zero" is not a good choice for numerics, so that was in a sense broken, too.
|
I agree, but the fix doesn't help a lot with the restrictions above.
Additionally I see the problem that you get different results with a MathIEEE program, depending on the version of the library (OS3.2 or not). For example the AmigaOS native version of vbcc is linked with -lmieee, so it might calculate constants differently on compile time under 3.2 than on run time under a different version.
But then it just adds to the initial problem of this thread, that also FPUs of cross-compilers seem to round differently...
Quote:
|
The only other option would be to set the FPCR on each call to the right precision. I can already hear the speed geeks cry.
|
Really no option, unfortunately.
Quote:
|
Concerning ANSI/C, isn't it so that computing (float) (op) (float) should be implemented as (round to float) ((double) (op) (double)), that is, your code should rather perform the computation in double precision, and then convert the result to single precision upon storing it in the variable?
|
Not sure. Have to check if ANSI/ISO-C allows that. Looking at the generated code I would say that (op) is done in the format of the arguments.
Quote:
|
Interesting, thanks for your findings. pow() is indeed a rather troublesome function (and for that reason also new in 3.9). Most folks just implement it by a^b=exp(b*log(a)), but this is not a good solution as it is numerically not stable. I don't recall what Hauser's code does here.
|
pow() is most likely not a problem of Hauser's Softfloat Library, as it has to be provided by support code. I'm using a modified version of e_pow.c from Sun Microsystem's portable mlib. It definitely does a more sophisticated approach than "a^b=exp(b*log(a))".
