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

fallback randn/randexp for AbstractFloat #44714

Merged
merged 6 commits into from
Apr 7, 2022
Merged

Conversation

stevengj
Copy link
Member

@stevengj stevengj commented Mar 23, 2022

This PR adds a fallback randn using the Marsaglia polar method (a Box–Muller variant) and a fallback randexp using -log(rand()) so that these two functions work for any AbstractFloat type implementing rand(rng, T).

As noted in #96 (comment), this gives us support for packages like DecFP.jl (closes JuliaMath/DecFP.jl#80) and DoubleFloats.jl for free.

Also closes #17629 for BigFloat (though a specialized method ala #35111 should be faster).

@andreasnoack, is there a way to test this using https://github.com/JuliaRandom/RNGTest.jl or similar?

@stevengj stevengj added the randomness Random number generation and the Random stdlib label Mar 23, 2022
@stevengj
Copy link
Member Author

stevengj commented Mar 23, 2022

(It's possible that the original Box–Muller algorithm is better here than the Marsaglia variant. It depends on the cost of a cos call compared to a 22% probability of 2 additional rand calls?)

I did some benchmarks with Float64, Dec64 (from DecFP.jl), and BigFloat, and the Marsaglia algorithm always seems to be faster (on average) than Box–Muller. The exception to the rule is Dec128, where Marsaglia is 20% slower for some reason.

@stevengj
Copy link
Member Author

The SparseArrays test checks that randn(BigFloat) throws a MethodError: https://github.com/JuliaSparse/SparseArrays.jl/blob/97189f3f36cce0a44210853bccfa4163e30678f3/test/sparse.jl#L3001

How can I update this test given that this is in another repo?

@DilumAluthge
Copy link
Member

Hmmm. Maybe first make a PR to SparseArrays.jl to remove that test? And then bump SparseArrays.jl on Julia master.

# this method is necessary to prevent rand(rng::AbstractRNG, X) from
# recursively constructing nested Sampler types.
Sampler(T::Type{<:AbstractRNG}, sp::Sampler, r::Repetition) =
throw(MethodError(Sampler, (T, sp, r)))
Copy link
Member Author

@stevengj stevengj Mar 24, 2022

Choose a reason for hiding this comment

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

Note that this is a minor change, but a MethodError seems to be what is really desired here (and is what we normally throw in methods we've only defined to prevent infinite dispatch loops).

(Otherwise, some downstream tests in SparseArrays fail because randn now falls through to rand and throws a different error, so this seems like the least-breaking change.)

@vtjnash
Copy link
Member

vtjnash commented Mar 28, 2022

need to address CI error:

Error in testset SparseArrays/sparse:
Test Failed at C:\buildbot\worker-tabularasa\tester_win64\build\share\julia\stdlib\v1.9\SparseArrays\test\sparse.jl:2976
  Expression: sprandn(T, 5, 5, 0.5)
    Expected: MethodError
  No exception thrown

@DilumAluthge
Copy link
Member

Needs #44731

@stevengj
Copy link
Member Author

Rebased, so it should include the fix #44731 now.

@stevengj stevengj added merge me PR is reviewed. Merge when all tests are passing triage This should be discussed on a triage call and removed merge me PR is reviewed. Merge when all tests are passing labels Apr 4, 2022
@stevengj
Copy link
Member Author

stevengj commented Apr 4, 2022

Okay to merge?

@stevengj
Copy link
Member Author

stevengj commented Apr 7, 2022

Unrelated CI failures.

@oscardssmith oscardssmith merged commit 244ada3 into master Apr 7, 2022
@oscardssmith oscardssmith deleted the fallback_randn branch April 7, 2022 23:31
@oscardssmith oscardssmith removed the triage This should be discussed on a triage call label May 26, 2022
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
randomness Random number generation and the Random stdlib
Projects
None yet
Development

Successfully merging this pull request may close these issues.

Lacking randn() Consider adding BigFloat support for randn/randexp
6 participants