Add support for local RNG seed in Haar measure matrix generation and reproducibility tests

src/Ginibre.jl

@@ -1,5 +1,6 @@
 export rand, Ginibre
 import Base.rand
+import Random: AbstractRNG, GLOBAL_RNG
 
 """
     Ginibre(β::Int, N::Int) <: ContinuousMatrixDistribution
@@ -30,30 +31,38 @@
 end
 
 """
-    rand(W::Ginibre)
+    rand(rng::AbstractRNG, W::Ginibre)
 
 Samples a matrix from the Ginibre ensemble.
 
 For `β = 1,2,4`, generates matrices randomly sampled from the real, complex, and quaternion
 Ginibre ensemble, respectively.
 """
-function rand(W::Ginibre)
+function rand(rng::AbstractRNG, W::Ginibre)
     beta, n = W.beta, W.N
     if beta==1
-        randn(n,n)
+        randn(rng,n,n)
     elseif beta==2
-        randn(n,n)+im*randn(n,n)
+        randn(rng,n,n)+im*randn(rng,n,n)
     elseif beta==4
-        Q0=randn(n,n)
-        Q1=randn(n,n)
-        Q2=randn(n,n)
-        Q3=randn(n,n)
+        Q0=randn(rng,n,n)
+        Q1=randn(rng,n,n)
+        Q2=randn(rng,n,n)
+        Q3=randn(rng,n,n)
         [Q0+im*Q1 Q2+im*Q3;-Q2+im*Q3 Q0-im*Q1]
     else 
         error(string("beta = ", beta, " not implemented"))
     end
 end
 
+
+rand(W::Ginibre)=rand(GLOBAL_RNG, W)
+
+
+
+
+
+
 function jpdf(Z::AbstractMatrix{z}) where {z<:Complex}
     pi^(size(Z,1)^2)*exp(-trace(Z'*Z))
 end

src/HaarMeasure.jl

@@ -1,6 +1,6 @@
 #TODO implement O(n^2) method
 """
-    rand(W::Haar, n::Int)
+    rand(rng::AbstractRNG, W::Haar, n::Int)
 
 Computes samples of real or complex Haar matrices of size `n`×`n`.
 
@@ -24,19 +24,19 @@
 - Edelman and Rao, 2005
 - Mezzadri, 2006, math-ph/0609050
 """
-function rand(W::Haar, n::Int, doCorrection::Int=1)
+function rand(rng::AbstractRNG, W::Haar, n::Int, doCorrection::Int=1)
     beta = W.beta
-    M=rand(Ginibre(beta,n))
+    M=rand(rng,Ginibre(beta,n))
     q,r=qr(M)
     if doCorrection==0
         q
     elseif doCorrection==1
         if beta==1
-            L = sign.(rand(n).-0.5)
+            L = sign.(rand(rng,n).-0.5)
         elseif beta==2
-            L = exp.(im*rand(n)*2pi)
+            L = exp.(im*rand(rng,n)*2pi)
         elseif beta==4
-            L = exp.(im*rand(2n)*2pi)
+            L = exp.(im*rand(rng,2n)*2pi)
         else
             error(string("beta = ",beta, " not implemented."))
         end
@@ -54,6 +54,10 @@
     end
 end
 
+rand(W::Haar, n::Int, doCorrection::Int=1)= rand(GLOBAL_RNG,W,n,doCorrection)
+
+
+
 #A utility method to check if the piecewise correction is needed
 #This checks the R part of the QR factorization; if correctly done,
 #the diagonals are all chi variables so are non-negative

test/Haar.jl

@@ -1,6 +1,7 @@
 using RandomMatrices
 using LinearAlgebra: I, tr, Diagonal, QRPackedQ
 using Test
+using Random
 
 @testset "Haar" begin
 
@@ -29,4 +30,23 @@ for T in (Float64, ComplexF64)
     @test C*A' ≈ C*A2'
 end
 
+
+@testset "Local RNG reproducibility" begin
+    rng1 = MersenneTwister(1234)
+    rng2 = MersenneTwister(1234)
+
+    A = rand(rng1, Haar(1), N)
+    B = rand(rng2, Haar(1), N)
+
+    @test A ≈ B
+
+    # Confirm that rand still works without passing an RNG
+    C = rand(Haar(1), N)
+    @test size(C) == (N, N)
+end
+
+
 end # testset
+
+
+