Best practices to ensure determinstic behaviour when random numbers are used

General Usage
1

Hi!
I have been playing around with some number theory algorithms like Rabin-Shallit’s algorithm (see LagrangeFourSquares.jl) for finding a four-square representation (i.e., given n find a,b,c,d such that n=a^2+b^2+c^2+d^2). However, this tasks has two issues:

  • There could be many solutions. For example, 28=1^2+1^2+1^2+5^2=1^2+3^2+3^2+3^2.
  • The algorithm that I am working with uses a random number generator to obtain some intermediary guesses.

So my concern is that this function is not deterministic. I am wondering if other Julians had a similar issue and how they dealt with it. Right now, I basically initialize an RNG with the same seed whenever the function is called to keep it deterministic. Is that what is usually done?

The closest analogy I can think of is finding a divisor of n by randomly testing numbers below \sqrt{n}:

using Random

f(n)=begin
    rng=MersenneTwister(1234)
    isqrtn=isqrt(n)
    while true 
        m=rand(rng,2:isqrtn)
        mod(n,m)==0 && return m
    end
end

This appears to be deterministic but if the RNG was updated, it would break tests like f(10)==2.

2

StableRNGs.jl may be useful

3

See the discussion of exactly this topic in the manual: Random Numbers: Reproducibility

If you really want a deterministic algorithm, of course, then one could question whether you should be using pseudo-random numbers at all, as opposed to a low-discrepancy sequence (for a Quasi Monte-Carlo algorithm), e.g. using Sobol.jl.