Joe wrote:
BTW, AFAICT, DJGPP does not have a Log1P function, so I wasn't able to verify the above anyway. I didn't know the values to use for the coefficients "Log1PCoefs[]", so I couldn't try the suggested Log1P code, either. Although
??? You don't have to guess the coefficients, they are initialized in the code. (If you don't want to rebuild the RTS just to test it (I didn't either), move the "to begin do" stuff into your main program, and import the GPC unit (because of the SplitReal function).)
Emil presented a power series which I assume represents these coefficients, I don't know whether these need to be "tweaked". If I understand this correctly, it is an infinite series. However, AFAIK, when you apply a finite version of such a series, you need to "tweak" it's values to get good accuracy (no, I don't know how to do this).
I remember experimenting years ago with the power series for "sine" or "cosine" and being amazed at how a small finite "tweaked" power series was much more accurate than the "conventional" power series, no matter how many terms I calculated for it (after a while, the more terms you calculate, the greater the error you introduce). I presume the "log1p" power series also has this effect(?)
Not sure what sort of "tweaking" you have in mind. Anyway, my Log1P implementation tries to use as few terms of the power series as possible, and it sums them in reverse order. I don't know if this is optimal, but it appeared to work with acceptable precision: the largest error I found was 2 ulp's when compiled without optimization, and 1 ulp when compiled with -O2, both on an i586.
Emil