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

widen signature of randn!(::MersenneTwister, ::Array{Float64}) ? #49388

Open
bjarthur opened this issue Apr 17, 2023 · 5 comments
Open

widen signature of randn!(::MersenneTwister, ::Array{Float64}) ? #49388

bjarthur opened this issue Apr 17, 2023 · 5 comments
Labels
randomness Random number generation and the Random stdlib

Comments

@bjarthur
Copy link
Contributor

why should this performance optimization not include all AbstractArrays{Float64}s ? i haven't tested any changes to the Random stdlib, but simply looking at the code i think it should work for AbstractArrays too, like this: randfun!(rng::MersenneTwister, A::AbstractArray{Float64}).

the reason i ask, is that for testing purposes i need to generate the same random numbers with and without physical units using Unitful.jl. to do so, one currently must strip off the units using ustrip, which uses reinterpret, and hence returns a ReinterpretArray, which then dispatches to the non-optimized code, leading to different random numbers.

julia> using Unitful, Random

julia> len=15;

julia> rng = MersenneTwister();

julia> a = zeros(len);

julia> Random.seed!(rng, 1);

julia> randn!(rng, ustrip(a))
15-element Vector{Float64}:
  0.2972879845354616
  0.3823959677906078
 -0.5976344767282311
 -0.01044524463737564
 -0.839026854388764
  0.31111133849833383
  2.2950878238373105
  0.40839583832475224
  0.2290095549097807
 -2.2670863488005306
  0.5299655761667461
  0.43142152642291204
  0.5837082875687786
  0.9632716050381906
  0.45879095505371686

julia> ua = zeros(len) * 1u"s";

julia> Random.seed!(rng, 1);

julia> randn!(rng, ustrip(ua))
15-element reinterpret(Float64, ::Vector{Quantity{Float64, 𝐓, Unitful.FreeUnits{(s,), 𝐓, nothing}}}):
  0.2972879845354616
  0.3823959677906078
 -0.5976344767282311
 -0.01044524463737564
 -0.839026854388764
  0.31111133849833383
  2.2950878238373105
 -2.2670863488005306
  0.5299655761667461
  0.43142152642291204
  0.5837082875687786
  0.9632716050381906
  0.45879095505371686
 -0.5223367574215084
  0.40839583832475224

cc: @rfourquet

@Seelengrab
Copy link
Contributor

Using @inbounds on a type you don't own is generally not safe to do, so limiting this method to the (Base-owned) Array is the correct way to have the use be safe.

@bjarthur
Copy link
Contributor Author

the fallback method, which inputs AbstractArray (as well as AbstractRNG), uses @inbounds too, so widening the input signature of the optimized method for Float64 (and MersenneTwister) to include all AbstractArray not just Array would take no particular risk.

@oscardssmith
Copy link
Member

that's bad. We should fix that.

@Seelengrab
Copy link
Contributor

Seelengrab commented Apr 18, 2023

The risk is still there - another method having that same problem doesn't make it any more safe to do.

@udohjeremiah udohjeremiah added the randomness Random number generation and the Random stdlib label Apr 19, 2023
@rfourquet
Copy link
Member

rfourquet commented May 17, 2023

I personally think it's better to use @inbounds even for non-Array arrays, but the reason it was limited to Array here is that the effects were not measured. If no optimized method exists for rand!(floatarray), it's not clear that randn!(floatarray) would benefit from the optimization of the linked PR. So someone would need to make some experiments.

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

No branches or pull requests

5 participants