Skip to content

Conversation

@oscardssmith
Copy link
Member

Not sure if this is a regression, but it's definitely bad. Also, fixes a problems I found in expm1(::Float32) with very small inputs.

@oscardssmith oscardssmith added maths Mathematical functions bugfix This change fixes an existing bug embarrassing-bugfix Whoops! float16 labels May 19, 2021
@oscardssmith oscardssmith added this to the 1.7 milestone May 19, 2021
@jarlebring
Copy link
Contributor

A bit off-topic: Why is this file mixing so much between using evalpoly and exthorner? Isn't horner to prefer over standard poly evaluation for higher degrees? Also, isn't it better to use poly evaluations with less multiplications, such as Paterson-Stockmayer? Or newer poly evaluations methods such as the one by Sastre: https://doi.org/10.1016/j.laa.2017.11.010

@oscardssmith
Copy link
Member Author

evalpoly uses horner's method. exthorner is horner with extended precision by using error free transforms. The link you gave is for evaluating matrix polynomials which is very different from scalar polynomials. There are theoretically more efficient methods for polynomials than horner's method (estirn's for one), but they are only worth-while for polynomials of relatively high degree (12 or more roughly). For small polynomials, you can't do anything more efficient than the series of fma instructions that evalpoly produces.

@jarlebring
Copy link
Contributor

jarlebring commented May 19, 2021

Ok. Thanks.

You have a point that the paper concerns matrix polynomials. I agree that it is more complicated at this low-level, but I thought a rule-of-thumb of computation cost for scalar is the same: Multiplication is much more expensive than addition, and division (which is not present in the paper I referenced) is even more expensive.

For small polynomials, you can't do anything more efficient than the series of fma instructions that evalpoly produces.

Okay. I can accept this directly assuming "small polynomials" means the degree is max (say) 4 or 5. The file has polys with degrees 8 (or maybe higher). Horner eval of poly of degree 8 requires 8 multiplications. With 8 multiplications we can obtain a poly of degree 25 using Paterson-Stockmayer, see figure 1 in http://eprints.maths.manchester.ac.uk/2676/1/fasi18.pdf

I don't understand why a Paterson-Stockmayer approach would be less favorable for fma's.

The estrin scheme also seems to require quite a few multiplications.

@jarlebring
Copy link
Contributor

jarlebring commented May 19, 2021

Ping @mfasi who is an expert on Paterson-Stockmayer. What do you think about the use of PS for scalar polynomial evaluation?

@jarlebring
Copy link
Contributor

jarlebring commented May 19, 2021

Okay. On second thought maybe my comment is less useful for this than I first thought. You have to count also coefficient * argument multiplication as a regular multiplication and the PS is not as competitive.

@jarlebring
Copy link
Contributor

It seems it is possible to get less fma instructions with a PS approach for a polynomial of degree 8. I am happy to discuss it further in a better place: https://discourse.julialang.org/t/fmas-and-poly-eval/61421

@oscardssmith
Copy link
Member Author

Can I get review on this? This fixes a fairly major bug (expm(small_x)), and adds expm(::Float16) which definitely should exist.

@oscardssmith oscardssmith requested review from simonbyrne and removed request for musm May 24, 2021 05:52
Copy link
Contributor

@musm musm left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In general, looks fine to me. I did not verify coefficients, but I trust they're good presumably generated via same method (Remez.jl).

@oscardssmith
Copy link
Member Author

What's buildkite? I haven't seen it before.

Co-authored-by: Mustafa M <mus-m@outlook.com>
@oscardssmith oscardssmith removed the request for review from simonbyrne May 24, 2021 18:53
@oscardssmith
Copy link
Member Author

Merging in a day or 2 sans objections.

@oscardssmith oscardssmith merged commit 2d1f16e into JuliaLang:master May 26, 2021
@oscardssmith oscardssmith deleted the float16-exp branch May 26, 2021 20:10
@oscardssmith oscardssmith removed this from the 1.7 milestone May 27, 2021
shirodkara pushed a commit to shirodkara/julia that referenced this pull request Jun 9, 2021
add expm1(::Float16) and tests for expm1(-floatmax(T))

Co-authored-by: Mustafa M <mus-m@outlook.com>
Co-authored-by: Kristoffer Carlsson <kcarlsson89@gmail.com>
johanmon pushed a commit to johanmon/julia that referenced this pull request Jul 5, 2021
add expm1(::Float16) and tests for expm1(-floatmax(T))

Co-authored-by: Mustafa M <mus-m@outlook.com>
Co-authored-by: Kristoffer Carlsson <kcarlsson89@gmail.com>
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

bugfix This change fixes an existing bug embarrassing-bugfix Whoops! float16 maths Mathematical functions

Projects

None yet

Development

Successfully merging this pull request may close these issues.

4 participants