Improve besselj accuracy for large arguments for orders -1, 0, 1

[jwscook]

Related to #37, which gives examples for order 11 and argument 1e9.

It looks as though the large argument expansions for besselj0 and besselj1 are also inaccurate. I made a gist comparing the differences with SpecialFunctions.jl over a large range of orders and magnitude of argument.

The gist produces this figure indicating that the errors blow up for orders -1, 0 and 1 and arguments > 1e16.
BesseljComparison.pdf

[heltonmc]

Very much thank you for this report!!!! And the details are so helpful.

@oscardssmith this looks related to cos_sum failing for large arguments. When you have time could you have a look at that function when x>1e16.

Bessels.jl/src/misc.jl

Lines 4 to 20 in 3bd0109

	function cos_sum(x, xn)
	s = x + xn
	n, r = reduce_pi02_med(s)
	lo = r.lo - ((s - x) - xn)
	hi = r.hi + lo
	y = Base.Math.DoubleFloat64(hi, r.hi-hi+lo)
	n = n&3
	if n == 0
	return Base.Math.cos_kernel(y)
	elseif n == 1
	return -Base.Math.sin_kernel(y)
	elseif n == 2
	return -Base.Math.cos_kernel(y)
	else
	return Base.Math.sin_kernel(y)
	end
	end

On a quick inspection..

julia> cos_sum(1e20, 1.2)
-3.5964366282318264e27

julia> cos(1e20 + 1.2)
0.7639704044417283

[heltonmc]

#13 (comment)

Relevant discussion on where this was introduced. Sure enough as Oscar said as long as x<1e16. I’m sorry for this I should have had tests for very large inputs.

[jwscook]

Great! Such fast work :)

[heltonmc]

Thank you very much for the bug report! It is much appreciated. Hopefully this is fixed on master and will tag a new release soon. If the situation persists please let me know and can reopen!