aboutsummaryrefslogtreecommitdiffstats
path: root/lib/stdlib/src/rand.erl
diff options
context:
space:
mode:
authorRaimo Niskanen <[email protected]>2018-06-18 10:49:41 +0200
committerRaimo Niskanen <[email protected]>2018-09-13 13:35:14 +0200
commitae074a6d4ed9a056d6687b2736be28a9dcc0649d (patch)
tree5d5ef8c97460c73ab8a6d7eb7a61e524051eb61f /lib/stdlib/src/rand.erl
parentc4408318c97be75915b3dff23a933bbc4b5c2046 (diff)
downloadotp-ae074a6d4ed9a056d6687b2736be28a9dcc0649d.tar.gz
otp-ae074a6d4ed9a056d6687b2736be28a9dcc0649d.tar.bz2
otp-ae074a6d4ed9a056d6687b2736be28a9dcc0649d.zip
Improve seeding methods
Diffstat (limited to 'lib/stdlib/src/rand.erl')
-rw-r--r--lib/stdlib/src/rand.erl298
1 files changed, 178 insertions, 120 deletions
diff --git a/lib/stdlib/src/rand.erl b/lib/stdlib/src/rand.erl
index 3a32f40fcd..fdf9709633 100644
--- a/lib/stdlib/src/rand.erl
+++ b/lib/stdlib/src/rand.erl
@@ -133,9 +133,10 @@
exrop | exs1024s | exro928ss | exsp | exs64 | exsplus | exs1024.
-type alg() :: builtin_alg() | atom().
-type export_state() :: {alg(), alg_state()}.
+-type seed() :: [integer()] | integer() | {integer(), integer(), integer()}.
-export_type(
[builtin_alg/0, alg/0, alg_handler/0, alg_state/0,
- state/0, export_state/0]).
+ state/0, export_state/0, seed/0]).
-export_type(
[exrop_state/0, exs1024_state/0, exro928_state/0, exsplus_state/0,
exs64_state/0]).
@@ -236,12 +237,12 @@ export_seed() ->
end.
-spec export_seed_s(State :: state()) -> export_state().
-export_seed_s({#{type:=Alg}, Seed}) -> {Alg, Seed}.
+export_seed_s({#{type:=Alg}, AlgState}) -> {Alg, AlgState}.
%% seed(Alg) seeds RNG with runtime dependent values
%% and return the NEW state
-%% seed({Alg,Seed}) setup RNG with a previously exported seed
+%% seed({Alg,AlgState}) setup RNG with a previously exported seed
%% and return the NEW state
-spec seed(
@@ -253,11 +254,11 @@ seed(Alg) ->
-spec seed_s(
AlgOrStateOrExpState :: builtin_alg() | state() | export_state()) ->
state().
-seed_s({AlgHandler, _Seed} = State) when is_map(AlgHandler) ->
+seed_s({AlgHandler, _AlgState} = State) when is_map(AlgHandler) ->
State;
-seed_s({Alg0, Seed}) ->
- {Alg,_SeedFun} = mk_alg(Alg0),
- {Alg, Seed};
+seed_s({Alg, AlgState}) when is_atom(Alg) ->
+ {AlgHandler,_SeedFun} = mk_alg(Alg),
+ {AlgHandler,AlgState};
seed_s(Alg) ->
seed_s(Alg, {erlang:phash2([{node(),self()}]),
erlang:system_time(),
@@ -266,21 +267,15 @@ seed_s(Alg) ->
%% seed/2: seeds RNG with the algorithm and given values
%% and returns the NEW state.
--spec seed(
- Alg :: builtin_alg(),
- Seed :: {integer(), integer(), integer()}) ->
- state().
-seed(Alg0, S0) ->
- seed_put(seed_s(Alg0, S0)).
+-spec seed(Alg :: builtin_alg(), Seed :: seed()) -> state().
+seed(Alg, Seed) ->
+ seed_put(seed_s(Alg, Seed)).
--spec seed_s(
- Alg :: builtin_alg(),
- Seed :: {integer(), integer(), integer()}) ->
- state().
-seed_s(Alg0, S0) ->
- {Alg, Seed} = mk_alg(Alg0),
- AS = Seed(S0),
- {Alg, AS}.
+-spec seed_s(Alg :: builtin_alg(), Seed :: seed()) -> state().
+seed_s(Alg, Seed) ->
+ {AlgHandler,SeedFun} = mk_alg(Alg),
+ AlgState = SeedFun(Seed),
+ {AlgHandler,AlgState}.
%%% uniform/0, uniform/1, uniform_s/1, uniform_s/2 are all
%%% uniformly distributed random numbers.
@@ -290,8 +285,8 @@ seed_s(Alg0, S0) ->
-spec uniform() -> X :: float().
uniform() ->
- {X, Seed} = uniform_s(seed_get()),
- _ = seed_put(Seed),
+ {X, State} = uniform_s(seed_get()),
+ _ = seed_put(State),
X.
%% uniform/1: given an integer N >= 1,
@@ -300,8 +295,8 @@ uniform() ->
-spec uniform(N :: pos_integer()) -> X :: pos_integer().
uniform(N) ->
- {X, Seed} = uniform_s(N, seed_get()),
- _ = seed_put(Seed),
+ {X, State} = uniform_s(N, seed_get()),
+ _ = seed_put(State),
X.
%% uniform_s/1: given a state, uniform_s/1
@@ -650,6 +645,14 @@ mk_alg(exro928ss) ->
-opaque exs64_state() :: uint64().
+exs64_seed(L) when is_list(L) ->
+ [R] = seed64_nz(1, L),
+ R;
+exs64_seed(A) when is_integer(A) ->
+ [R] = seed64(1, ?MASK(64, A)),
+ R;
+%%
+%% Traditional integer triplet seed
exs64_seed({A1, A2, A3}) ->
{V1, _} = exs64_next((?MASK(32, A1) * 4294967197 + 1)),
{V2, _} = exs64_next((?MASK(32, A2) * 4294967231 + 1)),
@@ -676,6 +679,14 @@ exs64_next(R) ->
-dialyzer({no_improper_lists, exsplus_seed/1}).
+exsplus_seed(L) when is_list(L) ->
+ [S0,S1] = seed58_nz(2, L),
+ [S0|S1];
+exsplus_seed(X) when is_integer(X) ->
+ [S0,S1] = seed58(2, ?MASK(64, X)),
+ [S0|S1];
+%%
+%% Traditional integer triplet seed
exsplus_seed({A1, A2, A3}) ->
{_, R1} = exsplus_next(
[?MASK(58, (A1 * 4294967197) + 1)|
@@ -751,6 +762,12 @@ exsplus_jump(S, [AS0|AS1], J, N) ->
-opaque exs1024_state() :: {list(uint64()), list(uint64())}.
+exs1024_seed(L) when is_list(L) ->
+ {seed64_nz(16, L), []};
+exs1024_seed(X) when is_integer(X) ->
+ {seed64(16, ?MASK(64, X)), []};
+%%
+%% Seed from traditional triple, remain backwards compatible
exs1024_seed({A1, A2, A3}) ->
B1 = ?MASK(21, (?MASK(21, A1) + 1) * 2097131),
B2 = ?MASK(21, (?MASK(21, A2) + 1) * 2097133),
@@ -931,32 +948,17 @@ exs1024_jump({L, RL}, AS, JL, J, N, TN) ->
-opaque exro928_state() :: {list(uint58()), list(uint58())}.
-
-%%% %% Seed raw words
-%%% exro928_seed({S0,S1,S2,S3,S4,S5,S6,S7,S8,S9,S10,S11,S12,S13,S14,S15})
-%%% when S0 bor S1 bor S2 bor S3 bor S4 bor S5 bor S6 bor S7 bor
-%%% S8 bor S9 bor S10 bor S11 bor S12 bor S13 bor S14 bor S15 > 0,
-%%% S0 bor S1 bor S2 bor S3 bor S4 bor S5 bor S6 bor S7 bor
-%%% S8 bor S9 bor S10 bor S11 bor S12 bor S13 bor S14 bor S15 < 1 bsl 58 ->
-%%% {[S0,S1,S2,S3,S4,S5,S6,S7,S8,S9,S10,S11,S12,S13,S14,S15], []};
-%%% %%
-%%% %% Seed from one 64-bit integer through splitmix
-%%% exro928_seed(X) when is_integer(X), 0 =< X, X =< 1 bsl 64 ->
-%%% {exro928_seed(X, 16),[]};
+exro928_seed(L) when is_list(L) ->
+ {seed58_nz(16, L), []};
+exro928_seed(X) when is_integer(X) ->
+ {seed58(16, ?MASK(64, X)), []};
%%
-%% Seed from traditional triple - splitmix mixed with the 3 integers
+%% Seed from traditional integer triple - mix into splitmix
exro928_seed({A1, A2, A3}) ->
{S0, X0} = seed58(?MASK(64, A1)),
- {S1, X1} = seed58(?MASK(64, A2 bxor X0)),
- {S2, X2} = seed58(?MASK(64, A3 bxor X1)),
- {[S0,S1,S2|exro928_seed(X2, 13)],[]}.
-%%
-%% Splitmix seed the rest of the state words
-exro928_seed(_X, 0) ->
- [];
-exro928_seed(X, N) ->
- {S, NewX} = seed58(X),
- [S|exro928_seed(NewX, N-1)].
+ {S1, X1} = seed58(?MASK(64, A2) bxor X0),
+ {S2, X2} = seed58(?MASK(64, A3) bxor X1),
+ {[S0,S1,S2|seed58(13, X2)], []}.
%% Update the state and calculate output word
@@ -1038,65 +1040,6 @@ exro928_jump_2pow20(SR) ->
16#415FEFCC78FF45E, 16#72CF460728C2FAF]).
%% =====================================================================
-%% Polynomial jump with a jump constant word list,
-%% high bit in each word marking top of word,
-%% SR is a {Forward, Reverse} queue tuple with Forward never empty
-%% =====================================================================
-
-polyjump({Ss, Rs} = SR, NextState, JumpConst) ->
- %% Create new state accumulator T
- Ts = lists:duplicate(length(Ss) + length(Rs), 0),
- polyjump(SR, NextState, JumpConst, Ts).
-%%
-%% Foreach jump word
-polyjump(_SR, _NextState, [], Ts) ->
- %% Return new calculated state
- {Ts, []};
-polyjump(SR, NextState, [J|Js], Ts) ->
- polyjump(SR, NextState, Js, Ts, J).
-%%
-%% Foreach bit in jump word until top bit
-polyjump(SR, NextState, Js, Ts, 1) ->
-
- polyjump(SR, NextState, Js, Ts);
-polyjump({Ss, Rs} = SR, NextState, Js, Ts, J) when J =/= 0 ->
- NewSR = NextState(SR),
- NewJ = J bsr 1,
- case ?MASK(1, J) of
- 0 ->
- polyjump(NewSR, NextState, Js, Ts, NewJ);
- 1 ->
- %% Xor this state onto T
- polyjump(NewSR, NextState, Js, xorzip_sr(Ts, Ss, Rs), NewJ)
- end.
-
-xorzip_sr([], [], undefined) ->
- [];
-xorzip_sr(Ts, [], Rs) ->
- xorzip_sr(Ts, lists:reverse(Rs), undefined);
-xorzip_sr([T|Ts], [S|Ss], Rs) ->
- [T bxor S|xorzip_sr(Ts, Ss, Rs)].
-
-%% =====================================================================
-
-format_jumpconst58(String) ->
- ReOpts = [{newline,any},{capture,all_but_first,binary},global],
- {match,Matches} = re:run(String, "0x([a-zA-Z0-9]+)", ReOpts),
- format_jumcons58_matches(lists:reverse(Matches), 0).
-
-format_jumcons58_matches([], J) ->
- format_jumpconst58_value(J);
-format_jumcons58_matches([[Bin]|Matches], J) ->
- NewJ = (J bsl 64) bor binary_to_integer(Bin, 16),
- format_jumcons58_matches(Matches, NewJ).
-
-format_jumpconst58_value(0) ->
- ok;
-format_jumpconst58_value(J) ->
- io:format("16#~s,~n", [integer_to_list(?MASK(58, J) bor ?BIT(58), 16)]),
- format_jumpconst58_value(J bsr 58).
-
-%% =====================================================================
%% exrop PRNG: Xoroshiro116+
%%
%% Reference URL: http://xorshift.di.unimi.it/
@@ -1164,6 +1107,15 @@ format_jumpconst58_value(J) ->
-opaque exrop_state() :: nonempty_improper_list(uint58(), uint58()).
-dialyzer({no_improper_lists, exrop_seed/1}).
+
+exrop_seed(L) when is_list(L) ->
+ [S0,S1] = seed58_nz(2, L),
+ [S0|S1];
+exrop_seed(X) when is_integer(X) ->
+ [S0,S1] = seed58(2, ?MASK(64, X)),
+ [S0|S1];
+%%
+%% Traditional integer triplet seed
exrop_seed({A1, A2, A3}) ->
[_|S1] =
exrop_next_s(
@@ -1227,14 +1179,34 @@ exrop_jump([S__0|S__1] = _S, S0, S1, J, Js) ->
end.
%% =====================================================================
-%% 58-bit seeder; lowest 58 bits of SplitMix64, zeros skipped
-%%
-%% uint64_t splitmix64_next() {
-%% uint64_t z = (x += 0x9e3779b97f4a7c15);
-%% z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9;
-%% z = (z ^ (z >> 27)) * 0x94d049bb133111eb;
-%% return z ^ (z >> 31);
-%% }
+%% Mask and fill state list, ensure not all zeros
+%% =====================================================================
+
+seed58_nz(N, Ss) ->
+ seed_nz(N, Ss, 58, false).
+
+seed64_nz(N, Ss) ->
+ seed_nz(N, Ss, 64, false).
+
+seed_nz(_N, [], _M, false) ->
+ erlang:error(zero_seed);
+seed_nz(0, [_|_], _M, _NZ) ->
+ erlang:error(too_many_seed_integers);
+seed_nz(0, [], _M, _NZ) ->
+ [];
+seed_nz(N, [], M, true) ->
+ [0|seed_nz(N - 1, [], M, true)];
+seed_nz(N, [S|Ss], M, NZ) ->
+ if
+ is_integer(S) ->
+ R = ?MASK(M, S),
+ [R|seed_nz(N - 1, Ss, M, NZ orelse R =/= 0)];
+ true ->
+ erlang:error(non_integer_seed)
+ end.
+
+%% =====================================================================
+%% Splitmix seeders, lowest bits of SplitMix64, zeros skipped
%% =====================================================================
seed58(0, _X) ->
@@ -1244,16 +1216,102 @@ seed58(N, X) ->
[Z|seed58(N - 1, NewX)].
%%
seed58(X_0) ->
- X = ?MASK(64, X_0 + 16#9e3779b97f4a7c15),
- Z_0 = ?MASK(64, (X bxor (X bsr 30)) * 16#bf58476d1ce4e5b9),
- Z_1 = ?MASK(64, (Z_0 bxor (Z_0 bsr 27)) * 16#94d049bb133111eb),
- case ?MASK(58, Z_1 bxor (Z_1 bsr 31)) of
+ {Z0,X} = splitmix64_next(X_0),
+ case ?MASK(58, Z0) of
0 ->
seed58(X);
Z ->
- {Z, X}
+ {Z,X}
+ end.
+
+seed64(0, _X) ->
+ [];
+seed64(N, X) ->
+ {Z,NewX} = seed64(X),
+ [Z|seed64(N - 1, NewX)].
+%%
+seed64(X_0) ->
+ {Z,X} = ZX = splitmix64_next(X_0),
+ if
+ Z =:= 0 ->
+ seed64(X);
+ true ->
+ ZX
+ end.
+
+%% The SplitMix64 generator:
+%%
+%% uint64_t splitmix64_next() {
+%% uint64_t z = (x += 0x9e3779b97f4a7c15);
+%% z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9;
+%% z = (z ^ (z >> 27)) * 0x94d049bb133111eb;
+%% return z ^ (z >> 31);
+%% }
+%%
+splitmix64_next(X_0) ->
+ X = ?MASK(64, X_0 + 16#9e3779b97f4a7c15),
+ Z_0 = ?MASK(64, (X bxor (X bsr 30)) * 16#bf58476d1ce4e5b9),
+ Z_1 = ?MASK(64, (Z_0 bxor (Z_0 bsr 27)) * 16#94d049bb133111eb),
+ {?MASK(64, Z_1 bxor (Z_1 bsr 31)),X}.
+
+%% =====================================================================
+%% Polynomial jump with a jump constant word list,
+%% high bit in each word marking top of word,
+%% SR is a {Forward, Reverse} queue tuple with Forward never empty
+%% =====================================================================
+
+polyjump({Ss, Rs} = SR, NextState, JumpConst) ->
+ %% Create new state accumulator T
+ Ts = lists:duplicate(length(Ss) + length(Rs), 0),
+ polyjump(SR, NextState, JumpConst, Ts).
+%%
+%% Foreach jump word
+polyjump(_SR, _NextState, [], Ts) ->
+ %% Return new calculated state
+ {Ts, []};
+polyjump(SR, NextState, [J|Js], Ts) ->
+ polyjump(SR, NextState, Js, Ts, J).
+%%
+%% Foreach bit in jump word until top bit
+polyjump(SR, NextState, Js, Ts, 1) ->
+ polyjump(SR, NextState, Js, Ts);
+polyjump({Ss, Rs} = SR, NextState, Js, Ts, J) when J =/= 0 ->
+ NewSR = NextState(SR),
+ NewJ = J bsr 1,
+ case ?MASK(1, J) of
+ 0 ->
+ polyjump(NewSR, NextState, Js, Ts, NewJ);
+ 1 ->
+ %% Xor this state onto T
+ polyjump(NewSR, NextState, Js, xorzip_sr(Ts, Ss, Rs), NewJ)
end.
+xorzip_sr([], [], undefined) ->
+ [];
+xorzip_sr(Ts, [], Rs) ->
+ xorzip_sr(Ts, lists:reverse(Rs), undefined);
+xorzip_sr([T|Ts], [S|Ss], Rs) ->
+ [T bxor S|xorzip_sr(Ts, Ss, Rs)].
+
+%% =====================================================================
+
+format_jumpconst58(String) ->
+ ReOpts = [{newline,any},{capture,all_but_first,binary},global],
+ {match,Matches} = re:run(String, "0x([a-zA-Z0-9]+)", ReOpts),
+ format_jumcons58_matches(lists:reverse(Matches), 0).
+
+format_jumcons58_matches([], J) ->
+ format_jumpconst58_value(J);
+format_jumcons58_matches([[Bin]|Matches], J) ->
+ NewJ = (J bsl 64) bor binary_to_integer(Bin, 16),
+ format_jumcons58_matches(Matches, NewJ).
+
+format_jumpconst58_value(0) ->
+ ok;
+format_jumpconst58_value(J) ->
+ io:format("16#~s,~n", [integer_to_list(?MASK(58, J) bor ?BIT(58), 16)]),
+ format_jumpconst58_value(J bsr 58).
+
%% =====================================================================
%% Ziggurat cont
%% =====================================================================