File 2952-Remove-inferior-35-bit-PRNGs.patch of Package erlang
From 0c8130f3cb573daf2d98d46b50b3f732e2dfa267 Mon Sep 17 00:00:00 2001
From: Raimo Niskanen <raimo@erlang.org>
Date: Mon, 2 May 2022 13:56:56 +0200
Subject: [PATCH 2/6] Remove inferior 35-bit PRNGs
---
lib/stdlib/doc/src/rand.xml | 151 -------------------------
lib/stdlib/src/rand.erl | 106 +-----------------
lib/stdlib/test/rand_SUITE.erl | 196 ++-------------------------------
3 files changed, 10 insertions(+), 443 deletions(-)
diff --git a/lib/stdlib/doc/src/rand.xml b/lib/stdlib/doc/src/rand.xml
index 884251ab02..22135d3443 100644
--- a/lib/stdlib/doc/src/rand.xml
+++ b/lib/stdlib/doc/src/rand.xml
@@ -408,14 +408,6 @@ tests. We suggest to use a sign test to extract a random Boolean value.</pre>
<name name="uint64"/>
<desc><p>0 .. (2^64 - 1)</p></desc>
</datatype>
- <datatype>
- <name name="mcg35_state"/>
- <desc><p>1 .. ((2^35 - 31) - 1)</p></desc>
- </datatype>
- <datatype>
- <name name="lcg35_state"/>
- <desc><p>0 .. (2^35 - 1)</p></desc>
- </datatype>
<datatype>
<name name="mwc59_state"/>
<desc><p>1 .. ((16#1ffb072 * 2^29 - 1) - 1)</p></desc>
@@ -893,149 +885,6 @@ end.</pre>
</p>
</desc>
</func>
- <func>
- <name name="mcg35" arity="1" since="OTP 25.0"/>
- <fsummary>Return a new generator state / number.</fsummary>
- <desc>
- <p>
- Returns a generated Pseudo Random Number which is also
- the new generated state, <c><anno>X1</anno></c>,
- according to a classical Multiplicative Congruential Generator
- (a.k.a Multiplicative Linear Congruential Generator,
- Lehmer random number generator,
- Park-Miller random number generator).
- </p>
- <p>
- This generator uses the modulus 2^35 - 31
- and the multiplication constant 185852
- from the paper "Tables of Linear Congruential Generators
- of different sizes and good lattice structure" by
- Pierre L'Ecuyer (1997) and they are selected
- for performance to keep the computation under
- the Erlang bignum limit.
- </p>
- <p>
- The generator may be written as
- <c><anno>X1</anno> = (185852*<anno>X0</anno>) rem ((1 bsl 35)-31)</c>,
- but the properties of the chosen constants has allowed
- an optimization of the otherwise expensive <c>rem</c> operation.
- </p>
- <p>
- On a typical 64 bit Erlang VM this generator executes
- in just below 10% (1/10) of the time
- for the default algorithm in this module.
- </p>
- <note>
- <p>
- This generator is only suitable for insensitive
- special niche applications since it has
- a short period (2^35 - 32),
- few bits (under 35),
- is not a power of 2 generator (range 1 .. (2^35 - 32)),
- offers no help in generating numbers on a specified range,
- and so on.
- </p>
- <p>
- But for pseudo random load distribution
- and such it might be useful, since it is very fast.
- It normally beats even the known-to-be-fast trick
- <c>erlang:phash2(erlang:unique_integer())</c>.
- </p>
- </note>
- </desc>
- </func>
- <func>
- <name name="mcg35_seed" arity="0" since="OTP 25.0"/>
- <name name="mcg35_seed" arity="1" since="OTP 25.0"/>
- <fsummary>Create a generator state.</fsummary>
- <desc>
- <p>
- Returns a generator state <c><anno>X</anno></c>.
- It is set it to <c><anno>S</anno></c>, after folding it
- into the state's range, if out of range.
- </p>
- <p>
- Without <c><anno>S</anno></c>,
- the generator state is created as for
- <seemfa marker="#seed_s/1"><c>seed_s(atom())</c></seemfa>.
- </p>
- </desc>
- </func>
- <func>
- <name name="lcg35" arity="1" since="OTP 25.0"/>
- <fsummary>Return a new generator state / number.</fsummary>
- <desc>
- <p>
- Returns a generated Pseudo Random Number which is also
- the new generated state, <c><anno>X1</anno></c>,
- according to a classical Linear Congruential Generator,
- a power of 2 mixed congruential generator.
- </p>
- <p>
- This generator uses the modulus 2^35
- and the multiplication constant 15319397
- from the paper "Tables of Linear Congruential Generators
- of different sizes and good lattice structure" by
- Pierre L'Ecuyer (1997) and they are selected
- for performance to keep the computation under
- the Erlang bignum limit.
- The addition constant has been selected to
- 15366142135 (it has to be odd) which looks
- more interesting than simply 1.
- </p>
- <p>
- The generator may be written as
- <c><anno>X1</anno> = ((15319397*<anno>X0</anno>) + 15366142135) band ((1 bsl 35)-1)</c>.
- </p>
- <p>
- On a typical 64 bit Erlang VM this generator executes
- in just below 7% (1/15) of the time
- for the default algorithm in this module,
- which is the fastest generator the author has seen.
- It can hardly be beaten by even a BIF
- implementation of any generator since the execution
- time is close to the overhead of a BIF call.
- </p>
- <note>
- <p>
- This generator is only suitable for insensitive
- special niche applications since it has
- a short period (2^35), few bits (35),
- offers no help in generating numbers on a specified range,
- and has, among others, the known statistical artifact
- that the lowest bit simply alternates,
- the next to lowest has a period of 4,
- and so on, so it is only the highest bits
- that achieve any form of statistical quality.
- </p>
- <p>
- But for pseudo random load distribution
- and such it might be useful, since it is extremely fast.
- The <seemfa marker="#mcg35/1"><c>mcg35/1</c></seemfa>
- generator above has less statistical artifacts,
- but instead it has other peculiarities since it is not
- a power of 2 generator.
- </p>
- </note>
- </desc>
- </func>
- <func>
- <name name="lcg35_seed" arity="0" since="OTP 25.0"/>
- <name name="lcg35_seed" arity="1" since="OTP 25.0"/>
- <fsummary>Create a generator state.</fsummary>
- <desc>
- <p>
- Returns a generator state <c><anno>X</anno></c>.
- It is set it to <c><anno>S</anno></c>, after folding it
- into the state's range, if out of range.
- </p>
- <p>
- Without <c><anno>S</anno></c>,
- the generator state is created as for
- <seemfa marker="#seed_s/1"><c>seed_s(atom())</c></seemfa>.
- </p>
- </desc>
- </func>
<func>
<name name="mwc59" arity="1" since="OTP 25.0"/>
<fsummary>Return a new generator state.</fsummary>
diff --git a/lib/stdlib/src/rand.erl b/lib/stdlib/src/rand.erl
index aa228d6eaa..217b1d14b5 100644
--- a/lib/stdlib/src/rand.erl
+++ b/lib/stdlib/src/rand.erl
@@ -38,8 +38,6 @@
%% Utilities
-export([exsp_next/1, exsp_jump/1, splitmix64_next/1,
- mcg35/1, mcg35_seed/0, mcg35_seed/1,
- lcg35/1, lcg35_seed/0, lcg35_seed/1,
mwc59/1, mwc59_value/1, mwc59_full_value/1, mwc59_float/1,
mwc59_seed/0, mwc59_seed/1]).
@@ -151,7 +149,7 @@
[exsplus_state/0, exro928_state/0, exrop_state/0, exs1024_state/0,
exs64_state/0, mwc59_state/0, dummy_state/0]).
-export_type(
- [uint58/0, uint64/0, splitmix64_state/0, mcg35_state/0, lcg35_state/0]).
+ [uint58/0, uint64/0, splitmix64_state/0]).
%% =====================================================================
%% Range macro and helper
@@ -1450,108 +1448,6 @@ dummy_seed({A1, A2, A3}) ->
?MASK(58, Z3).
-%% =====================================================================
-%% mcg35 PRNG: Multiplicative Linear Congruential Generator
-%%
-%% Parameters by Pierre L'Ecuyer from the paper:
-%% TABLES OF LINEAR CONGRUENTIAL GENERATORS
-%% OF DIFFERENT SIZES AND GOOD LATTICE STRUCTURE
-%%
-%% X1 = (A * X0) rem M
-%%
-%% A = 185852
-%% M = 2^B - D
-%% B = 35, D = 31
-%%
-%% X cannot, and never will be 0.
-%% Take X1 as output value with range 1..(M-1).
-%%
-%% Use this generator for speed, not for quality.
-%% The calculation should avoid bignum operations,
-%% so A * X0 should not become a bignum.
-%% M and A has been selected accordingly.
-%% The period is M - 1 since 0 cannot occur.
-%%
-%% =====================================================================
--define(MCG35_A, (185852)).
--define(MCG35_B, (35)).
--define(MCG35_D, (31)).
--define(MCG35_M, (?BIT(?MCG35_B) - ?MCG35_D)).
-
--type mcg35_state() :: 1..(?MCG35_M-1).
-
--spec mcg35(X0 :: mcg35_state()) -> X1 :: mcg35_state().
-mcg35(X0) -> % when is_integer(X0), 1 =< X0, X0 < ?MCG35_M ->
- %% The mask operation on the input tricks the JIT into
- %% realizing that all following operations does not
- %% need bignum handling. The suggested guard test above
- %% could have had the same effect but it did not, and, alas,
- %% needs much more native code to execute than a 'band'.
- X = ?MCG35_A * ?MASK(?MCG35_B, X0),
- %% rem M = rem (2^B - D), optimization to avoid 'rem'
- X1 = ?MASK(?MCG35_B, X) + ?MCG35_D*(X bsr ?MCG35_B),
- if
- X1 < ?MCG35_M -> X1;
- true -> X1 - ?MCG35_M
- end.
-
--spec mcg35_seed() -> X :: mcg35_state().
-mcg35_seed() ->
- {A1, A2, A3} = default_seed(),
- mcg35_seed(seed64_int([A1, A2, A3])).
-
--spec mcg35_seed(S :: integer()) -> X :: mcg35_state().
-mcg35_seed(S) when is_integer(S) ->
- mod(?MCG35_M - 1, S) + 1.
-
-%% =====================================================================
-%% lcg35 PRNG: Multiplicative Linear Congruential Generator
-%%
-%% Parameters by Pierre L'Ecuyer from the paper:
-%% TABLES OF LINEAR CONGRUENTIAL GENERATORS
-%% OF DIFFERENT SIZES AND GOOD LATTICE STRUCTURE
-%%
-%% X1 = (A * X0 + C) rem M
-%%
-%% M = 2^35, A = 15319397, C = 1
-%%
-%% Choosing C = 15366142135, which is an odd value
-%% about 2^35 / sqrt(5) gives a perhaps nicer value in
-%% the sequence after 0 (-> C).
-%%
-%% Take X1 as output value with range 1..(M-1).
-%%
-%% Use this generator for speed, not for quality.
-%% The calculation should avoid bignum operations,
-%% so A * X0 + C should not become a bignum.
-%% M, A and C has been selected accordingly.
-%% The period is M.
-%%
-%% =====================================================================
--define(LCG35_A, (15319397)).
--define(LCG35_C, (15366142135)). % (1 bsl 35) / sqrt(5), odd
--define(LCG35_B, 35). % Number of bits
-
--type lcg35_state() :: 0..?MASK(?LCG35_B).
-
--spec lcg35(X0 :: lcg35_state()) -> X1 :: lcg35_state().
-lcg35(X0) -> % when is_integer(X0), 0 =< X0, X0 =< ?MASK(?LCG35_B) ->
- %% The mask operation on the input tricks the JIT into
- %% realizing that all following operations does not
- %% need bignum handling. The suggested guard test above
- %% could have had the same effect but it did not, and, alas,
- %% needs much more native code to execute than a 'band'.
- ?MASK(?LCG35_B, ?LCG35_A * ?MASK(?LCG35_B, X0) + ?LCG35_C).
-
--spec lcg35_seed() -> X :: lcg35_state().
-lcg35_seed() ->
- {A1, A2, A3} = default_seed(),
- lcg35_seed(seed64_int([A1, A2, A3])).
-
--spec lcg35_seed(S :: integer()) -> X :: lcg35_state().
-lcg35_seed(S) when is_integer(S) ->
- ?MASK(?LCG35_B, S).
-
%% =====================================================================
%% mcg58 PRNG: Multiply With Carry generator
%%
diff --git a/lib/stdlib/test/rand_SUITE.erl b/lib/stdlib/test/rand_SUITE.erl
index 93d7b2914a..54363b2761 100644
--- a/lib/stdlib/test/rand_SUITE.erl
+++ b/lib/stdlib/test/rand_SUITE.erl
@@ -35,8 +35,9 @@ all() ->
[seed, interval_int, interval_float,
bytes_count,
api_eq,
- mcg35_api, mcg35_rem, lcg35_api, mwc59_api,
- exsp_next_api, exsp_jump_api, splitmix64_next_api,
+ mwc59_api,
+ exsp_next_api, exsp_jump_api,
+ splitmix64_next_api,
reference,
{group, basic_stats},
{group, distr_stats},
@@ -207,53 +208,6 @@ api_eq_1(S00) ->
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
--define(MCG35_M, ((1 bsl 35) - 31)).
-
-%% Verify mcg35 behaviour
-%%
-mcg35_api(Config) when is_list(Config) ->
- mcg35_api(1, 1000000).
-
-mcg35_api(X0, 0) ->
- X = 30203473714,
- {X, X} = {X0, X},
- ok;
-mcg35_api(X0, N)
- when is_integer(X0), 1 =< X0, X0 < ?MCG35_M ->
- mcg35_api(rand:mcg35(X0), N - 1).
-
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
-%% Verify the 'rem' optimization in mcg35
-%%
-mcg35_rem(Config) when is_list(Config) ->
- {_, X0} = rand:seed_s(dummy),
- mcg35_rem((X0 rem (?MCG35_M - 1)) + 1, 10000000).
-
-mcg35_rem(_X0, 0) ->
- ok;
-mcg35_rem(X0, N) ->
- X1 = (185852 * X0) rem ?MCG35_M,
- {X1, X1, X0, N} = {rand:mcg35(X0), X1, X0, N},
- mcg35_rem(X1, N - 1).
-
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
-%% Verify lcg35 behaviour
-%%
-lcg35_api(Config) when is_list(Config) ->
- lcg35_api(0, 1000000).
-
-lcg35_api(X0, 0) ->
- X = 19243690048,
- {X, X} = {X0, X},
- ok;
-lcg35_api(X0, N)
- when is_integer(X0), 0 =< X0, X0 < 1 bsl 35 ->
- lcg35_api(rand:lcg35(X0), N - 1).
-
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
%% Verify mwc59 behaviour
%%
mwc59_api(Config) when is_list(Config) ->
@@ -1153,35 +1107,6 @@ do_measure(Iterations) ->
Generator(Range, St0), Range, X, St1)
end
end, Algs, Iterations),
- _ =
- measure_1(
- fun (_Mod, _State) ->
- Range = 10000,
- fun (St0) ->
- St1 = rand:mcg35(St0),
- %% Just a 'rem' with slightly skewed distribution
- case (St1 rem Range) + 1 of
- R when is_integer(R), 0 < R, R =< Range ->
- St1
- end
- end
- end, {mcg35,mod}, Iterations,
- TMarkUniformRange10000, OverheadUniformRange1000),
- _ =
- measure_1(
- fun (_Mod, _State) ->
- Range = 10000,
- fun (St0) ->
- St1 = rand:lcg35(St0),
- %% Just a 'rem' with slightly skewed distribution
- case (St1 rem Range) + 1 of
- R when is_integer(R), 0 < R, R =< Range ->
- St1
- end
- end
- end,
- {lcg35,mod}, Iterations,
- TMarkUniformRange10000, OverheadUniformRange1000),
_ =
measure_1(
fun (_Mod, _State) ->
@@ -1286,20 +1211,6 @@ do_measure(Iterations) ->
end
end,
Algs, Iterations),
- _ =
- measure_1(
- fun (_Mod, _State) ->
- Range = 1 bsl 32,
- fun (St0) ->
- St1 = rand:lcg35(St0),
- case St1 bsr (35-32) of
- R when is_integer(R), 0 =< R, R < Range ->
- St1
- end
- end
- end,
- {lcg35,shift}, Iterations,
- TMarkUniform32Bit, OverheadUniform32Bit),
_ =
measure_1(
fun (_Mod, _State) ->
@@ -1440,32 +1351,6 @@ do_measure(Iterations) ->
Generator(Range, St0), Range, X, St1)
end
end, Algs ++ [dummy], Iterations),
- _ =
- measure_1(
- fun (_Mod, _State) ->
- Range = (1 bsl 35) - 31,
- fun (St0) ->
- case rand:mcg35(St0) of
- R when is_integer(R), 1 =< R, R =< Range ->
- R
- end
- end
- end,
- {mcg35,raw}, Iterations,
- TMarkUniformFullRange, OverheadUniformFullRange),
- _ =
- measure_1(
- fun (_Mod, _State) ->
- Range = 1 bsl 35,
- fun (St0) ->
- case rand:lcg35(St0) of
- R when is_integer(R), 0 =< R, R < Range ->
- R
- end
- end
- end,
- {lcg35,raw}, Iterations,
- TMarkUniformFullRange, OverheadUniformFullRange),
_ =
measure_1(
fun (_Mod, _State) ->
@@ -1570,18 +1455,18 @@ do_measure(Iterations) ->
_ =
measure_1(
fun (_Mod, _State) ->
- Range = 1 bsl 35,
+ Range = (16#7f17555 bsl 32) - 1,
fun (St0) ->
case
- put(lcg35_procdict,
- rand:lcg35(get(lcg35_procdict)))
+ put(mwc59_procdict,
+ rand:mwc59(get(mwc59_procdict)))
of
X when is_integer(X), 0 =< X, X < Range ->
St0
end
end
end,
- {lcg35,procdict}, Iterations,
+ {mwc59,procdict}, Iterations,
TMarkUniformFullRange, OverheadUniformFullRange),
%%
ct:pal("~nRNG uniform integer full range + 1 performance~n",[]),
@@ -1663,15 +1548,6 @@ do_measure(Iterations) ->
_ ->
Algs
end, Iterations),
- _ =
- measure_1(
- fun (_Mod, _State) ->
- fun (St0) ->
- ?CHECK_BYTE_SIZE(
- lcg35_bytes(ByteSize, St0), ByteSize, Bin, St1)
- end
- end, {lcg35,bytes}, Iterations,
- TMarkBytes1, OverheadBytes1),
_ =
measure_1(
fun (_Mod, _State) ->
@@ -1699,15 +1575,6 @@ do_measure(Iterations) ->
_ ->
Algs
end, Iterations div 50),
- _ =
- measure_1(
- fun (_Mod, _State) ->
- fun (St0) ->
- ?CHECK_BYTE_SIZE(
- lcg35_bytes(ByteSize2, St0), ByteSize2, Bin, St1)
- end
- end, {lcg35,bytes}, Iterations div 50,
- TMarkBytes2, OverheadBytes2),
_ =
measure_1(
fun (_Mod, _State) ->
@@ -1728,34 +1595,6 @@ do_measure(Iterations) ->
end
end,
Algs, Iterations),
- _ =
- measure_1(
- fun (_Mod, _State) ->
- fun (St0) ->
- St1 = rand:mcg35(St0),
- %% Too few bits for full mantissa
- case (St1 - 1) * (1/((1 bsl 35) - 31)) of
- R when is_float(R), 0.0 =< R, R < 1.0 ->
- St1
- end
- end
- end,
- {mcg35,float35}, Iterations,
- TMarkUniformFloat, OverheadUniformFloat),
- _ =
- measure_1(
- fun (_Mod, _State) ->
- fun (St0) ->
- St1 = rand:lcg35(St0),
- %% Too few bits for full mantissa
- case St1 * (1/(1 bsl 35)) of
- R when is_float(R), 0.0 =< R, R < 1.0 ->
- St1
- end
- end
- end,
- {lcg35,float35}, Iterations,
- TMarkUniformFloat, OverheadUniformFloat),
_ =
measure_1(
fun (_Mod, _State) ->
@@ -1916,13 +1755,9 @@ measure_init(Alg) ->
{rand, rand:seed(exsss)};
{Name, Tag} ->
case Name of
- lcg35 when Tag =:= procdict ->
- _ = put(lcg35_procdict, rand:lcg35_seed()),
+ mwc59 when Tag =:= procdict ->
+ _ = put(mwc59_procdict, rand:mwc59_seed()),
{rand, undefined};
- lcg35 ->
- {rand, rand:lcg35_seed()};
- mcg35 ->
- {rand, rand:mcg35_seed()};
mwc59 ->
{rand, rand:mwc59_seed()};
exsp ->
@@ -1948,19 +1783,6 @@ bytes_s(N, undefined = St) ->
%% crypto_bytes
{crypto:strong_rand_bytes(N), St}.
-lcg35_bytes(N, R0) ->
- lcg35_bytes(N, R0, <<>>).
-%%
-lcg35_bytes(N, R0, Bin) when N < 4 ->
- R1 = rand:lcg35(R0),
- Bits = N bsl 3,
- V = R1 bsr (35 - Bits),
- {<<Bin/binary, V:Bits>>, R1};
-lcg35_bytes(N, R0, Bin) ->
- R1 = rand:lcg35(R0),
- V = R1 bsr (35 - 32),
- lcg35_bytes(N - 4, R1, <<Bin/binary, V:32>>).
-
mwc59_bytes(N, R0) ->
mwc59_bytes(N, R0, <<>>).
%%
--
2.35.3