1 /*
2 * Licensed to the Apache Software Foundation (ASF) under one or more
3 * contributor license agreements. See the NOTICE file distributed with
4 * this work for additional information regarding copyright ownership.
5 * The ASF licenses this file to You under the Apache License, Version 2.0
6 * (the "License"); you may not use this file except in compliance with
7 * the License. You may obtain a copy of the License at
8 *
9 * http://www.apache.org/licenses/LICENSE-2.0
10 *
11 * Unless required by applicable law or agreed to in writing, software
12 * distributed under the License is distributed on an "AS IS" BASIS,
13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 * See the License for the specific language governing permissions and
15 * limitations under the License.
16 */
17
18 package org.apache.commons.rng.simple.internal;
19
20 import org.apache.commons.rng.UniformRandomProvider;
21 import org.apache.commons.rng.core.util.NumberFactory;
22
23 /**
24 * Utility for creating seeds.
25 */
26 final class SeedUtils {
27 /** The modulus {@code 256 % 9}. */
28 private static final int MOD_09 = 256 % 9;
29 /** The modulus {@code 256 % 10}. */
30 private static final int MOD_10 = 256 % 10;
31 /** The modulus {@code 256 % 11}. */
32 private static final int MOD_11 = 256 % 11;
33 /** The modulus {@code 256 % 12}. */
34 private static final int MOD_12 = 256 % 12;
35 /** The modulus {@code 256 % 13}. */
36 private static final int MOD_13 = 256 % 13;
37 /** The modulus {@code 256 % 14}. */
38 private static final int MOD_14 = 256 % 14;
39 /** The modulus {@code 256 % 15}. */
40 private static final int MOD_15 = 256 % 15;
41 /** The 16 hex digits in an array. */
42 private static final byte[] HEX_DIGIT_ARRAY = {0xf, 0xe, 0xd, 0xc, 0xb, 0xa, 0x9, 0x8,
43 0x7, 0x6, 0x5, 0x4, 0x3, 0x2, 0x1, 0x0};
44
45 /**
46 * Provider of unsigned 8-bit integers.
47 */
48 private static class UnsignedByteProvider {
49 /** Mask to extract the lowest 2 bits from an integer. */
50 private static final int MASK_2_BITS = 0x3;
51 /** Mask to extract the lowest 8 bits from an integer. */
52 private static final int MASK_8_BITS = 0xff;
53
54 /** Source of randomness. */
55 private final UniformRandomProvider rng;
56 /** The current 32-bits of randomness. */
57 private int bits;
58 /** The counter tracking the bits to extract. */
59 private int counter;
60
61 /**
62 * @param rng Source of randomness.
63 */
64 UnsignedByteProvider(UniformRandomProvider rng) {
65 this.rng = rng;
66 }
67
68 /**
69 * Return the next unsigned byte.
70 *
71 * @return Value in the range[0,255]
72 */
73 int nextUnsignedByte() {
74 // Every 4 bytes generate a new 32-bit value
75 final int count = counter & MASK_2_BITS;
76 counter++;
77 if (count == 0) {
78 bits = rng.nextInt();
79 // Consume the first byte
80 return bits & MASK_8_BITS;
81 }
82 // Consume a remaining byte (occurs 3 times)
83 bits >>>= 8;
84 return bits & MASK_8_BITS;
85 }
86 }
87
88 /**
89 * Class contains only static methods.
90 */
91 private SeedUtils() {}
92
93 /**
94 * Creates an {@code int} containing a permutation of 8 hex digits chosen from 16.
95 *
96 * @param rng Source of randomness.
97 * @return hex digit permutation.
98 */
99 static int createIntHexPermutation(UniformRandomProvider rng) {
100 final UnsignedByteProvider provider = new UnsignedByteProvider(rng);
101 return createUpperBitsHexPermutation(provider);
102 }
103
104 /**
105 * Creates a {@code long} containing a permutation of 8 hex digits chosen from 16 in
106 * the upper and lower 32-bits.
107 *
108 * @param rng Source of randomness.
109 * @return hex digit permutation.
110 */
111 static long createLongHexPermutation(UniformRandomProvider rng) {
112 final UnsignedByteProvider provider = new UnsignedByteProvider(rng);
113 // Extract upper bits and combine with a second sample
114 return NumberFactory.makeLong(createUpperBitsHexPermutation(provider),
115 createUpperBitsHexPermutation(provider));
116 }
117
118 /**
119 * Creates a {@code int} containing a permutation of 8 hex digits chosen from 16.
120 *
121 * @param provider Source of randomness.
122 * @return hex digit permutation.
123 */
124 private static int createUpperBitsHexPermutation(UnsignedByteProvider provider) {
125 // Compute a Fisher-Yates shuffle in-place on the 16 hex digits.
126 // Each digit is chosen uniformly from the remaining digits.
127 // The value is swapped with the current digit.
128
129 // The following is an optimised sampling algorithm that generates
130 // a uniform deviate in the range [0,n) from an unsigned byte.
131 // The expected number of samples is 256/n.
132 // This has a bias when n does not divide into 256. So samples must
133 // be rejected at a rate of (256 % n) / 256:
134 // n 256 % n Rejection rate
135 // 9 4 0.015625
136 // 10 6 0.0234375
137 // 11 3 0.01171875
138 // 12 4 0.015625
139 // 13 9 0.03515625
140 // 14 4 0.015625
141 // 15 1 0.00390625
142 // 16 0 0
143 // The probability of no rejections is 1 - (1-p15) * (1-p14) ... = 0.115
144
145 final byte[] digits = HEX_DIGIT_ARRAY.clone();
146
147 // The first digit has no bias. Get an index using a mask to avoid modulus.
148 int bits = copyToOutput(digits, 0, 15, provider.nextUnsignedByte() & 0xf);
149
150 // All remaining digits have a bias so use the rejection algorithm to find
151 // an appropriate uniform deviate in the range [0,n)
152 bits = copyToOutput(digits, bits, 14, nextUnsignedByteInRange(provider, MOD_15, 15));
153 bits = copyToOutput(digits, bits, 13, nextUnsignedByteInRange(provider, MOD_14, 14));
154 bits = copyToOutput(digits, bits, 12, nextUnsignedByteInRange(provider, MOD_13, 13));
155 bits = copyToOutput(digits, bits, 11, nextUnsignedByteInRange(provider, MOD_12, 12));
156 bits = copyToOutput(digits, bits, 10, nextUnsignedByteInRange(provider, MOD_11, 11));
157 bits = copyToOutput(digits, bits, 9, nextUnsignedByteInRange(provider, MOD_10, 10));
158 bits = copyToOutput(digits, bits, 8, nextUnsignedByteInRange(provider, MOD_09, 9));
159
160 return bits;
161 }
162
163
164 /**
165 * Get the next unsigned byte in the range {@code [0,n)} rejecting any over-represented
166 * sample value using the pre-computed modulus.
167 *
168 * <p>This algorithm is as per Lemire (2019) adapted for 8-bit arithmetic.</p>
169 *
170 * @param provider Provider of bytes.
171 * @param threshold Modulus threshold {@code 256 % n}.
172 * @param n Upper range (exclusive)
173 * @return Value.
174 * @see <a href="https://arxiv.org/abs/1805.10941">
175 * Lemire (2019): Fast Random Integer Generation in an Interval</a>
176 */
177 private static int nextUnsignedByteInRange(UnsignedByteProvider provider, int threshold, int n) {
178 // Rejection method using multiply by a fraction:
179 // n * [0, 2^8 - 1)
180 // ------------
181 // 2^8
182 // The result is mapped back to an integer and will be in the range [0, n)
183 int m;
184 do {
185 // Compute product of n * [0, 2^8 - 1)
186 m = n * provider.nextUnsignedByte();
187
188 // Test the sample uniformity.
189 } while ((m & 0xff) < threshold);
190 // Final divide by 2^8
191 return m >>> 8;
192 }
193
194 /**
195 * Copy the lower hex digit to the output bits. Swap the upper hex digit into
196 * the lower position. This is equivalent to a swap step of a Fisher-Yates shuffle on
197 * an array but the output of the shuffle are written to the bits.
198 *
199 * @param digits Digits.
200 * @param bits Bits.
201 * @param upper Upper index.
202 * @param lower Lower index.
203 * @return Updated bits.
204 */
205 private static int copyToOutput(byte[] digits, int bits, int upper, int lower) {
206 // Move the existing bits up and append the next hex digit.
207 // This is equivalent to swapping lower to upper.
208 final int newbits = bits << 4 | digits[lower] & 0xf;
209 // Swap upper to lower.
210 digits[lower] = digits[upper];
211 return newbits;
212 }
213 }