Re: [Mingw-users] expm1()

[Greg C. <chi...@co...>]

On 2005-01-19 04:01 AM, Danny Smith wrote:

> Greg Chicares wrote:
> 
>>log1p() would be much more useful if its companion expm1()
>>[C99 7.12.6.3] were provided as well. [...]
> 
> Off top of head, it would be quite easy to do
> 
> if (x < small_value)
> {
>   /* Taylor series expansion
>   return y such that additional term is smaller than  machine eps  */
> }
> else
>  return exp(x) - 1;
> 
> Not very fast,  but I'll bet I can find a good old Fortran IV routine that has
> been tested as accurate.

The hardware already provides F2XM1, which is both fast and accurate.

I'd provide a patch, except that it would take many months to get
assignment papers signed by my employer. But maybe I can do other
work that makes it trivial for someone else who doesn't have that
problem.

I see plenty of expm1() implementations online, but some don't meet
the FSF definition of 'free'; and those that do are GPL or LGPL,
which I believe MinGW won't use.

Maybe I could provide a narrative description of changes to
   mingw/runtime/mingwex/math/expl.c
that would be easy to translate. The intel instruction set actually
treats expm1() as the general case, and you have to undo the 'minus
one' part to implement plain exp(). It might be as simple as
eliminating the lines

        "fld1\n\t"               /* 4 1.0              */
        "faddp\n\t"		/* 3 2^(fract(x * log2(e))) */

but let me look into that.