[Ecls-list] log(complex) inaccurate for some values

Raymond Toy toy.raymond at gmail.com
Tue Apr 15 00:11:30 UTC 2008


Juan Jose Garcia-Ripoll wrote:
> On Mon, Apr 14, 2008 at 2:37 AM, Raymond Toy <toy.raymond at gmail.com> wrote:
>   
>>  Here's an outline for a one possible implementation of log1p.  Very likely
>> not as optimized as the C version, but should work with reasonable accuracy.
>>     
>
> Hi, looks nice. The version I found around, and which only activates
> when the C library does not support log1p (Windows, probably) is one I
> found googling around, and which compensates the errors in the
> ordinary logarithm.
>
> double log1p(double x)
> {
>         double u = 1.0 + x;
>         if (u == 1) {
>                 return 0.0;
>         } else {
>                 return (log(u) * x)/(u - 1.0);
>         }
> }
>
> But your implementation seems more stable. The only thing that looks
> expensive is the step for making "x" small. Wouldn't it suffice to use
> log1p(x) = series if |x| <= 1, and log1p(x) = log1p(1/x) + log(x)
> elsewhere?
>
>   
Hey, I like your algorithm better.  Clean and simple.  I don't 
understand how it works though, but that's ok.

The algorithm I gave was from Brent.  Making x small takes a little bit 
of time, but each step the number is halved, so it doesn't take too 
long.  But yeah, you can just use regular old log for x large enough.  I 
don't know where that breakpoint would be.

Ray





More information about the ecl-devel mailing list