This is the mail archive of the newlib@sourceware.org mailing list for the newlib project.


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]
Other format: [Raw text]

Re: Possible mistake in expf


On Thu, 11 Jun 2015, ethereal@ethv.net wrote:

> Hello everyone,
> 
> I believe I may have stumbled across a minor mistake in the expf()
> implementation in newlib's libm.
> 
> newlib/libm/mathfp/sf_pow.c:99 [1] is:
> 
> > x = exp (t);
> 
> However, I don't see why the function exp should be used instead of
> expf, as x itself is a single-precision float. At a quick glance I don't
> see any reason why this would benefit greatly from the increased
> precision the double-precision version of exp would provide, but I'm far
> from an expert on this sort of thing. If someone could confirm or refute
> my suspicions, that would be appreciated!

Since t is float, it won't actually help - but to get good results, t 
ought to be double (note that it's computed using the double version of 
log).  To get results within a few ulp for pow (x, y) as exp (y * log (x)), 
you need y * log (x) to have as many fractional bits as there are bits in 
the mantissa of the result type (that is, for results with large 
exponents, the mantissa precision needed for y * log (x) is about the sum 
of mantissa and exponent bits in the type in question, or about the number 
of bits in the representation for an IEEE type - so 32 bits intermediate 
precision is desirable for powf, rather than just the 24 bits of float).

-- 
Joseph S. Myers
joseph@codesourcery.com


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]