Code with apache-2 on /usr/src

[Adhemerval Zanella]

On 28/05/2018 19:18, Steve Kargl wrote:
> On Mon, May 28, 2018 at 06:12:13PM -0300, Adhemerval Zanella wrote:
>>
>>>> And is having a different algorithm for single and double prevision 
>>>> a blocker for a future patch proposal?
>>>
>>> No.  Given the comment in sinf.c that max ULP is 0.56072, I do note that
>>> the current implementation of sinf in lib/msun is more accurate (for
>>> interesting values of x).  I also looked at single/s_sincosf.c.  It is
>>> rather dubious to have 80+ digit numerical constants for a float, which
>>> at most has 9 relevant digits.
>>>
>>
>> Also keep in mind my initial idea is to propose patches only to expf, powf, 
>> logf, expf2, and log2f.
> 
> OK, so I peeked at expf.  Comment claims max ulp of 0.502.
> Exhaustive testing for normal numbers in relevent range for
> the current implementation of expf(x) shows
> 
> Interval tested: [-18,88.72]
> ULP: 0.90951,   x = -5.19804668e+00f, /* 0xc0a65666 */
> flt =  5.52735012e-03f, /* 0x3bb51ec6 */
> dbl =  5.5273505437686398e-03, /* 0x3f76a3d8, 0xdd1aae8e */
> 
> But, then one looks at implementation details.  msun's current
> implementation is written in terms of single precision; while
> the routine you're suggesting is written in terms of double_t.
> So, achieving 0.502 ULP is due to having 53-bits in intermediate
> results.  It appears that the algorithm of the suggested code 
> cannot easily be generalized to double and long double without
> implementing a multiple-precision routines.

This is indeed true for the default implementation, although the same repo
has alternative implementation that uses only float for expf, powf, and 
logf.  However, as far as I could evaluated, the optimized expf and powf
single version does not yield any gain over current FreeBSD version, only
for the logf I see some gains.

Do you see any issue about current approach of using intermediary double_t
for internal calculations?

> 
> Note, years ago, I submitted implementations for expf, exp, 
> ld80/expl, ld128/expl, logf, log, ld80/logl, and ld128/logl
> based on papers by PTP Tang [1,2].  My versions for single
> and double precision were not adopted even though these had
> better accuracy.  Either Bruce Evans improved or with Bruce's
> help I improved the ld80 and ld128 routines, which were added
> to msun.  I know Bruce fixed minor issues with the single 
> and double precision routines, but he has not submitted patches.
> 
> 1. PTP Tang, "Table-driven implementation of the exponential
>    function in IEEE floating-point arithmetic," ACM Trans. Math.
>    Soft., 15, 144-157 (1989).
> 
> 2. PTP Tang,  "Table-driven implementation of the logarithm
>    function in IEEE floating-point arithmetic," ACM Trans. Math.
>    Soft., 16, 378-400 (1990).
> 

Thanks for the links, do you recall why exactly your implementations were
not adopted? Do you think a similar proposal based on the arm repo would
be also rejected?