Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Improve performance for large arguments in besselj0 #76

Closed
heltonmc opened this issue Jan 30, 2023 · 1 comment
Closed

Improve performance for large arguments in besselj0 #76

heltonmc opened this issue Jan 30, 2023 · 1 comment

Comments

@heltonmc
Copy link
Member

The bug fix in #57 was fine for accuracy but this should be a lot better....

# fine performance

julia> @benchmark besselj0(x) setup=(x=130.0 + rand())
BenchmarkTools.Trial: 10000 samples with 999 evaluations.
 Range (min  max):  8.550 ns  34.869 ns  ┊ GC (min  max): 0.00%  0.00%
 Time  (median):     8.675 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   8.688 ns ±  0.572 ns  ┊ GC (mean ± σ):  0.00% ± 0.00%

                         ▄          █                         
  ▂▁▁▁▁▁▁▁▁▁▁▄▁▁▁▁▁▁▁▁▁▁▁█▁▁▁▁▁▁▁▁▁▁█▁▁▁▁▁▁▁▁▁▁▁▄▁▁▁▁▁▁▁▁▁▁▃ ▂
  8.55 ns        Histogram: frequency by time        8.76 ns <

 Memory estimate: 0 bytes, allocs estimate: 0.

## bad
julia> @benchmark besselj0(x) setup=(x=1e17 + rand())
BenchmarkTools.Trial: 10000 samples with 994 evaluations.
 Range (min  max):  30.600 ns  59.397 ns  ┊ GC (min  max): 0.00%  0.00%
 Time  (median):     30.768 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   30.845 ns ±  0.856 ns  ┊ GC (mean ± σ):  0.00% ± 0.00%

    ██                                                         
  ▃▇███▄▂▂▂▂▂▁▁▁▁▂▁▁▁▁▁▁▂▂▁▂▂▁▁▁▂▁▂▂▁▁▂▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▂▂▂▂▂ ▂
  30.6 ns         Histogram: frequency by time        33.5 ns <

 Memory estimate: 0 bytes, allocs estimate: 0.

Perhaps we could improve this whole function a little bit but I remember the problem being with rem_pio2 in Julia. Could revisit this.

At large arguments besselj0(x) ~ sqrt(2/pi) * cos(x - pi/4) / sqrt(x) so we should be able to do a lot better here then our current approach. We could potentially also improve the range of validity of this approximation by adding a term in the cos term besselj0(x) ~ sqrt(2/pi) * cos(x - pi/4 - 1/8x) / sqrt(x).
The polynomial outside of the cos quickly becomes one but the periodic nature of the cos function requires us to keep some more terms except for very large arguments.

No matter what though, computing cos(x + y) accurately is the important thing here. But we should be able to do much better because this is such an important function.

@heltonmc
Copy link
Member Author

Closed by #88

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant