Merge pull request #13815 from rfourquet/rand-bigrange-mergeversions

merge 2 versions (for GMP 5 & 6) of rand(::UnitRange{BigInt})

Author: rfourquet

base/gmp.jl

@@ -46,6 +46,7 @@ mutable struct BigInt <: Integer
     alloc::Cint
     size::Cint
     d::Ptr{Limb}
+
     function BigInt()
         b = new(zero(Cint), zero(Cint), C_NULL)
         MPZ.init!(b)
@@ -77,9 +78,10 @@ BigInt(x)
 function __init__()
     try
         if gmp_version().major != GMP_VERSION.major || gmp_bits_per_limb() != GMP_BITS_PER_LIMB
-            error(string("The dynamically loaded GMP library (version $(gmp_version()) with __gmp_bits_per_limb == $(gmp_bits_per_limb()))\n",
-                         "does not correspond to the compile time version (version $GMP_VERSION with __gmp_bits_per_limb == $GMP_BITS_PER_LIMB).\n",
-                         "Please rebuild Julia."))
+            msg = gmp_bits_per_limb() != GMP_BITS_PER_LIMB ? error : warn
+            msg(string("The dynamically loaded GMP library (version $(gmp_version()) with __gmp_bits_per_limb == $(gmp_bits_per_limb()))\n",
+                       "does not correspond to the compile time version (version $GMP_VERSION with __gmp_bits_per_limb == $GMP_BITS_PER_LIMB).\n",
+                       "Please rebuild Julia."))
         end
 
         ccall((:__gmp_set_memory_functions, :libgmp), Void,
@@ -114,6 +116,9 @@ gmpz(op::Symbol) = (Symbol(:__gmpz_, op), :libgmp)
 init!(x::BigInt) = (ccall((:__gmpz_init, :libgmp), Void, (mpz_t,), &x); x)
 init2!(x::BigInt, a) = (ccall((:__gmpz_init2, :libgmp), Void, (mpz_t, bitcnt_t), &x, a); x)
 
+realloc2!(x, a) = (ccall((:__gmpz_realloc2, :libgmp), Void, (mpz_t, bitcnt_t), &x, a); x)
+realloc2(a) = realloc2!(BigInt(), a)
+
 sizeinbase(a::BigInt, b) = Int(ccall((:__gmpz_sizeinbase, :libgmp), Csize_t, (mpz_t, Cint), &a, b))
 
 for op in (:add, :sub, :mul, :fdiv_q, :tdiv_q, :fdiv_r, :tdiv_r, :gcd, :lcm, :and, :ior, :xor)
@@ -196,6 +201,9 @@ cmp_si(a::BigInt, b) = ccall((:__gmpz_cmp_si, :libgmp), Cint, (mpz_t, Clong), &a
 cmp_ui(a::BigInt, b) = ccall((:__gmpz_cmp_ui, :libgmp), Cint, (mpz_t, Culong), &a, b) % Int
 cmp_d(a::BigInt, b) = ccall((:__gmpz_cmp_d, :libgmp), Cint, (mpz_t, Cdouble), &a, b) % Int
 
+mpn_cmp(a::Ptr{Limb}, b::Ptr{Limb}, c) = ccall((:__gmpn_cmp, :libgmp), Cint, (Ptr{Limb}, Ptr{Limb}, Clong), a, b, c)
+mpn_cmp(a::BigInt, b::BigInt, c) = mpn_cmp(a.d, b.d, c)
+
 get_str!(x, a, b::BigInt) = (ccall((:__gmpz_get_str,:libgmp), Ptr{Cchar}, (Ptr{Cchar}, Cint, mpz_t), x, a, &b); x)
 set_str!(x::BigInt, a, b) = ccall((:__gmpz_set_str, :libgmp), Cint, (mpz_t, Ptr{UInt8}, Cint), &x, a, b) % Int
 get_d(a::BigInt) = ccall((:__gmpz_get_d, :libgmp), Cdouble, (mpz_t,), &a)

base/random.jl

@@ -3,7 +3,7 @@
 module Random
 
 using Base.dSFMT
-using Base.GMP: GMP_VERSION, Limb, MPZ
+using Base.GMP: Limb, MPZ
 import Base: copymutable, copy, copy!, ==
 
 export srand,
@@ -592,23 +592,12 @@ for (T, U) in [(UInt8, UInt32), (UInt16, UInt32),
     end
 end
 
-if GMP_VERSION.major >= 6
-    struct RangeGeneratorBigInt <: RangeGenerator
-        a::BigInt             # first
-        m::BigInt             # range length - 1
-        nlimbs::Int           # number of limbs in generated BigInt's
-        mask::Limb            # applied to the highest limb
-    end
-
-else
-    struct RangeGeneratorBigInt <: RangeGenerator
-        a::BigInt             # first
-        m::BigInt             # range length - 1
-        limbs::Vector{Limb}   # buffer to be copied into generated BigInt's
-        mask::Limb            # applied to the highest limb
-
-        RangeGeneratorBigInt(a, m, nlimbs, mask) = new(a, m, Vector{Limb}(nlimbs), mask)
-    end
+struct RangeGeneratorBigInt <: RangeGenerator
+    a::BigInt         # first
+    m::BigInt         # range length - 1
+    nlimbs::Int       # number of limbs in generated BigInt's (z ∈ [0, m])
+    nlimbsmax::Int    # max number of limbs for z+a
+    mask::Limb        # applied to the highest limb
 end
 
 
@@ -619,7 +608,8 @@ function RangeGenerator(r::UnitRange{BigInt})
     nlimbs, highbits = divrem(nd, 8*sizeof(Limb))
     highbits > 0 && (nlimbs += 1)
     mask = highbits == 0 ? ~zero(Limb) : one(Limb)<<highbits - one(Limb)
-    return RangeGeneratorBigInt(first(r), m, nlimbs, mask)
+    nlimbsmax = max(nlimbs, abs(last(r).size), abs(first(r).size))
+    return RangeGeneratorBigInt(first(r), m, nlimbs, nlimbsmax, mask)
 end
 
 
@@ -649,30 +639,21 @@ function rand{T<:Integer, U<:Unsigned}(rng::AbstractRNG, g::RangeGeneratorInt{T,
     (unsigned(g.a) + rem_knuth(x, g.k)) % T
 end
 
-if GMP_VERSION.major >= 6
-    # mpz_limbs_write and mpz_limbs_finish are available only in GMP version 6
-    function rand(rng::AbstractRNG, g::RangeGeneratorBigInt)
-        x = BigInt()
-        while true
-            # note: on CRAY computers, the second argument may be of type Cint (48 bits) and not Clong
-            xd = MPZ.limbs_write!(x, g.nlimbs)
-            limbs = unsafe_wrap(Array, xd, g.nlimbs)
-            rand!(rng, limbs)
-            limbs[end] &= g.mask
-            MPZ.limbs_finish!(x, g.nlimbs)
-            x <= g.m && return MPZ.add!(x, g.a)
-        end
+function rand(rng::AbstractRNG, g::RangeGeneratorBigInt)
+    x = MPZ.realloc2(g.nlimbsmax*8*sizeof(Limb))
+    limbs = unsafe_wrap(Array, x.d, g.nlimbs)
+    while true
+        rand!(rng, limbs)
+        @inbounds limbs[end] &= g.mask
+        MPZ.mpn_cmp(x, g.m, g.nlimbs) <= 0 && break
     end
-else
-    function rand(rng::AbstractRNG, g::RangeGeneratorBigInt)
-        x = BigInt()
-        while true
-            rand!(rng, g.limbs)
-            g.limbs[end] &= g.mask
-            MPZ.import!(x, length(g.limbs), -1, sizeof(Limb), 0, 0, g.limbs)
-            x <= g.m && return MPZ.add!(x, g.a)
-        end
+    # adjust x.size (normally done by mpz_limbs_finish, in GMP version >= 6)
+    x.size = g.nlimbs
+    while x.size > 0
+        @inbounds limbs[x.size] != 0 && break
+        x.size -= 1
     end
+    MPZ.add!(x, g.a)
 end
 
 rand(rng::AbstractRNG, r::UnitRange{<:Union{Signed,Unsigned,BigInt,Bool}}) = rand(rng, RangeGenerator(r))