From a1ac90fe935b02ec66ccb39504ead5464a5383ab Mon Sep 17 00:00:00 2001 From: Dan Rose Date: Thu, 6 Oct 2022 17:48:34 -0500 Subject: [PATCH 1/3] random bigint --- src/bigints/random.nim | 34 ++++++++++++++++++++++++++++++++++ tests/trandom.nim | 22 ++++++++++++++++++++++ 2 files changed, 56 insertions(+) create mode 100644 src/bigints/random.nim create mode 100644 tests/trandom.nim diff --git a/src/bigints/random.nim b/src/bigints/random.nim new file mode 100644 index 0000000..b0b3936 --- /dev/null +++ b/src/bigints/random.nim @@ -0,0 +1,34 @@ +import ../bigints +import std/sequtils +import std/options +import std/random + +proc rand*(r: var Rand, x: Slice[BigInt]): BigInt = + assert(x.a <= x.b, "invalid range") + let + spread = x.b - x.a + # number of bits *not* including leading bit + nbits = spread.fastLog2 + # number of limbs to generate completely randomly + nFullLimbs = max(nbits div 32 - 1, 0) + # highest possible value of the top two limbs. + hi64Max = (spread shr (nFullLimbs*32)).toInt[:uint64].get() + while true: + # these limbs can be generated completely randomly + var limbs = newSeqWith(nFullLimbs, r.rand(uint32.low..uint32.high)) + # generate the top two limbs more carefully. This all but guarantees + # that the entire number is in the correct range + let hi64 = r.rand(uint64.low..hi64Max) + limbs.add(cast[uint32](hi64)) + limbs.add(cast[uint32](hi64 shr 32)) + result = initBigInt(limbs) + if result <= spread: + break + result += x.a + +proc rand*(r: var Rand, max: BigInt): BigInt = + rand(r, 0'bi..max) + +proc rand*(x: Slice[BigInt]): BigInt = rand(randState(), x) +proc rand*(max: BigInt): BigInt = rand(randState(), max) + diff --git a/tests/trandom.nim b/tests/trandom.nim new file mode 100644 index 0000000..b1998b4 --- /dev/null +++ b/tests/trandom.nim @@ -0,0 +1,22 @@ +import bigints +import bigints/random +import std/options + +block: # check uniformity + let lo = pow(10'bi, 90) + let hi = pow(10'bi, 100) + var total = 0'bi + let trials = 1000 + let nbuckets = 33 + var buckets = newSeq[int](nbuckets) + for x in 0..trials: + let r = rand(lo..hi) + doAssert(lo <= r) + doAssert(r <= hi) + total += r + let iBucket = (r-lo) div ((hi-lo) div initBigInt(nbuckets)) + buckets[iBucket.toInt[:int]().get()] += 1 + for x in buckets: + doAssert(trials/nbuckets*0.5 < float(x)) + doAssert(float(x) < trials/nbuckets*1.5) + From cef9ff542f58032e16693625bdb65976ef3c0613 Mon Sep 17 00:00:00 2001 From: Dan Rose Date: Fri, 7 Oct 2022 13:00:45 -0500 Subject: [PATCH 2/3] fixups doc comments change some procs to funcs add trandom to test task --- bigints.nimble | 2 +- src/bigints/random.nim | 12 ++++++++---- 2 files changed, 9 insertions(+), 5 deletions(-) diff --git a/bigints.nimble b/bigints.nimble index 5de9609..17b8d86 100644 --- a/bigints.nimble +++ b/bigints.nimble @@ -16,7 +16,7 @@ task test, "Test bigints": echo "testing " & backend & " backend" for gc in ["refc", "arc", "orc"]: echo " using " & gc & " GC" - for file in ["tbigints.nim", "tbugs.nim"]: + for file in ["tbigints.nim", "tbugs.nim", "trandom.nim"]: exec "nim r --hints:off --experimental:strictFuncs --backend:" & backend & " --gc:" & gc & " tests/" & file exec "nim doc --hints:off --backend:" & backend & " --gc:" & gc & " src/bigints.nim" diff --git a/src/bigints/random.nim b/src/bigints/random.nim index b0b3936..b63977c 100644 --- a/src/bigints/random.nim +++ b/src/bigints/random.nim @@ -3,7 +3,8 @@ import std/sequtils import std/options import std/random -proc rand*(r: var Rand, x: Slice[BigInt]): BigInt = +func rand*(r: var Rand, x: Slice[BigInt]): BigInt = + ## Return a random `BigInt`, within the given range, using the given state. assert(x.a <= x.b, "invalid range") let spread = x.b - x.a @@ -14,7 +15,7 @@ proc rand*(r: var Rand, x: Slice[BigInt]): BigInt = # highest possible value of the top two limbs. hi64Max = (spread shr (nFullLimbs*32)).toInt[:uint64].get() while true: - # these limbs can be generated completely randomly + # these limbs can be generated completely arbitrarily var limbs = newSeqWith(nFullLimbs, r.rand(uint32.low..uint32.high)) # generate the top two limbs more carefully. This all but guarantees # that the entire number is in the correct range @@ -26,9 +27,12 @@ proc rand*(r: var Rand, x: Slice[BigInt]): BigInt = break result += x.a -proc rand*(r: var Rand, max: BigInt): BigInt = +func rand*(r: var Rand, max: BigInt): BigInt = + ## Return a random non-negative `BigInt`, up to `max`, using the given state. rand(r, 0'bi..max) proc rand*(x: Slice[BigInt]): BigInt = rand(randState(), x) -proc rand*(max: BigInt): BigInt = rand(randState(), max) + ## Return a random `BigInt`, within the given range. +proc rand*(max: BigInt): BigInt = rand(randState(), max) + ## Return a random `BigInt`, up to `max`. From 537d0fedc9fe6e89fc5ead010aefa34fd74927b8 Mon Sep 17 00:00:00 2001 From: Dan Rose Date: Fri, 7 Oct 2022 13:25:40 -0500 Subject: [PATCH 3/3] fix for 1.4 compatibility --- src/bigints/random.nim | 7 ++++++- tests/trandom.nim | 6 +++--- 2 files changed, 9 insertions(+), 4 deletions(-) diff --git a/src/bigints/random.nim b/src/bigints/random.nim index b63977c..e4c7058 100644 --- a/src/bigints/random.nim +++ b/src/bigints/random.nim @@ -29,7 +29,12 @@ func rand*(r: var Rand, x: Slice[BigInt]): BigInt = func rand*(r: var Rand, max: BigInt): BigInt = ## Return a random non-negative `BigInt`, up to `max`, using the given state. - rand(r, 0'bi..max) + rand(r, 0.initBigInt..max) + +# backwards compatibility with 1.4 +when not defined(randState): + var state = initRand(777) + proc randState(): var Rand = state proc rand*(x: Slice[BigInt]): BigInt = rand(randState(), x) ## Return a random `BigInt`, within the given range. diff --git a/tests/trandom.nim b/tests/trandom.nim index b1998b4..c7f2b7e 100644 --- a/tests/trandom.nim +++ b/tests/trandom.nim @@ -3,9 +3,9 @@ import bigints/random import std/options block: # check uniformity - let lo = pow(10'bi, 90) - let hi = pow(10'bi, 100) - var total = 0'bi + let lo = pow(10.initBigInt, 90) + let hi = pow(10.initBigInt, 100) + var total = 0.initBigInt let trials = 1000 let nbuckets = 33 var buckets = newSeq[int](nbuckets)