This is an automated email from the ASF dual-hosted git repository. aherbert pushed a commit to branch master in repository https://gitbox.apache.org/repos/asf/commons-rng.git
commit 7314a5af2dabc45303fc602f69ab054c53c92e04 Author: Alex Herbert <aherb...@apache.org> AuthorDate: Mon Aug 5 11:00:41 2019 +0100 RNG-85: Added Middle-Square Weyl Sequence generator. --- .../core/source32/MiddleSquareWeylSequence.java | 122 +++++++++++++++++++++ .../org/apache/commons/rng/core/ProvidersList.java | 3 + .../source32/MiddleSquareWeylSequenceTest.java | 81 ++++++++++++++ 3 files changed, 206 insertions(+) diff --git a/commons-rng-core/src/main/java/org/apache/commons/rng/core/source32/MiddleSquareWeylSequence.java b/commons-rng-core/src/main/java/org/apache/commons/rng/core/source32/MiddleSquareWeylSequence.java new file mode 100644 index 0000000..b84cd47 --- /dev/null +++ b/commons-rng-core/src/main/java/org/apache/commons/rng/core/source32/MiddleSquareWeylSequence.java @@ -0,0 +1,122 @@ +/* + * Licensed to the Apache Software Foundation (ASF) under one or more + * contributor license agreements. See the NOTICE file distributed with + * this work for additional information regarding copyright ownership. + * The ASF licenses this file to You under the Apache License, Version 2.0 + * (the "License"); you may not use this file except in compliance with + * the License. You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0 + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + */ + +package org.apache.commons.rng.core.source32; + +import org.apache.commons.rng.core.util.NumberFactory; + +/** + * Middle Square Weyl Sequence Random Number Generator. + * + * <p>A fast all-purpose 32-bit generator. Memory footprint is 192 bits and the period is at least + * {@code 2^64}.</p> + * + * <p>Implementation is based on the paper + * <a href="https://arxiv.org/abs/1704.00358v3">Middle Square Weyl Sequence RNG</a>.</p> + * + * @see <a href="https://en.wikipedia.org/wiki/Middle-square_method">Middle Square Method</a> + * + * @since 1.3 + */ +public class MiddleSquareWeylSequence extends IntProvider { + /** Size of the seed array. */ + private static final int SEED_SIZE = 3; + + /** State of the generator. */ + private long x; + /** State of the Weyl sequence. */ + private long w; + /** + * Increment for the Weyl sequence. This must be odd to ensure a full period. + * + * <p>This is not final to support the restore functionality.</p> + */ + private long s; + + /** + * Creates a new instance. + * + * <p>Note: The generator output quality is highly dependent on the initial seed. + * If the generator is seeded poorly (for example with all zeros) it is possible the + * generator may output zero for many cycles before the internal state recovers randomness. + * The seed elements are used to set:</p> + * + * <ol> + * <li>The state of the generator + * <li>The state of the Weyl sequence + * <li>The increment of the Weyl sequence + * </ol> + * + * <p>The third element is set to odd to ensure a period of at least 2<sup>64</sup>. If the + * increment is of low complexity then the Weyl sequence does not contribute high quality + * randomness. It is recommended to use a permutation of 8 hex characters for the upper + * and lower 32-bits of the increment.</p> + * + * <p>The state of the generator is squared during each cycle. There is a possibility that + * different seeds can produce the same output, for example 0 and 2<sup>32</sup> produce + * the same square. This can be avoided by using the high complexity Weyl increment for the + * state seed element.</p> + * + * @param seed Initial seed. + * If the length is larger than 3, only the first 3 elements will + * be used; if smaller, the remaining elements will be automatically set. + */ + public MiddleSquareWeylSequence(long[] seed) { + if (seed.length < SEED_SIZE) { + final long[] tmp = new long[SEED_SIZE]; + fillState(tmp, seed); + setSeedInternal(tmp); + } else { + setSeedInternal(seed); + } + } + + /** + * Seeds the RNG. + * + * @param seed Seed. + */ + private void setSeedInternal(long[] seed) { + x = seed[0]; + w = seed[1]; + // Ensure the increment is odd to provide a maximal period Weyl sequence. + this.s = seed[2] | 1L; + } + + /** {@inheritDoc} */ + @Override + protected byte[] getStateInternal() { + return composeStateInternal(NumberFactory.makeByteArray(new long[] {x, w, s}), + super.getStateInternal()); + } + + /** {@inheritDoc} */ + @Override + protected void setStateInternal(byte[] state) { + final byte[][] c = splitStateInternal(state, SEED_SIZE * 8); + setSeedInternal(NumberFactory.makeLongArray(c[0])); + super.setStateInternal(c[1]); + } + + /** {@inheritDoc} */ + @Override + public int next() { + x *= x; + x += w += s; + return (int) (x = (x >>> 32) | (x << 32)); + } +} diff --git a/commons-rng-core/src/test/java/org/apache/commons/rng/core/ProvidersList.java b/commons-rng-core/src/test/java/org/apache/commons/rng/core/ProvidersList.java index 1c00bfb..9d4b6b5 100644 --- a/commons-rng-core/src/test/java/org/apache/commons/rng/core/ProvidersList.java +++ b/commons-rng-core/src/test/java/org/apache/commons/rng/core/ProvidersList.java @@ -34,6 +34,7 @@ import org.apache.commons.rng.core.source32.Well44497a; import org.apache.commons.rng.core.source32.Well44497b; import org.apache.commons.rng.core.source32.ISAACRandom; import org.apache.commons.rng.core.source32.MersenneTwister; +import org.apache.commons.rng.core.source32.MiddleSquareWeylSequence; import org.apache.commons.rng.core.source32.MultiplyWithCarry256; import org.apache.commons.rng.core.source32.KISSRandom; import org.apache.commons.rng.core.source32.PcgXshRr32; @@ -103,6 +104,8 @@ public final class ProvidersList { add(LIST32, new PcgXshRs32(new long[] {g.nextLong()})); add(LIST32, new PcgMcgXshRr32(g.nextLong())); add(LIST32, new PcgMcgXshRs32(g.nextLong())); + // Ensure a high complexity increment is used for the Weyl sequence + add(LIST32, new MiddleSquareWeylSequence(new long[] {g.nextLong(), g.nextLong(), 0xb5ad4eceda1ce2a9L})); // ... add more here. // "long"-based RNGs. diff --git a/commons-rng-core/src/test/java/org/apache/commons/rng/core/source32/MiddleSquareWeylSequenceTest.java b/commons-rng-core/src/test/java/org/apache/commons/rng/core/source32/MiddleSquareWeylSequenceTest.java new file mode 100644 index 0000000..7a6c590 --- /dev/null +++ b/commons-rng-core/src/test/java/org/apache/commons/rng/core/source32/MiddleSquareWeylSequenceTest.java @@ -0,0 +1,81 @@ +/* + * Licensed to the Apache Software Foundation (ASF) under one or more + * contributor license agreements. See the NOTICE file distributed with + * this work for additional information regarding copyright ownership. + * The ASF licenses this file to You under the Apache License, Version 2.0 + * (the "License"); you may not use this file except in compliance with + * the License. You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0 + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + */ +package org.apache.commons.rng.core.source32; + +import org.apache.commons.rng.core.RandomAssert; +import org.junit.Test; + +public class MiddleSquareWeylSequenceTest { + @Test + public void testReferenceCode() { + /* + * The data was generated using the following program based on the author's C code: + * https://mswsrng.wixsite.com/rand + * + * #include <stdio.h> + * #include <stdint.h> + * + * uint64_t x, w, s; + * + * inline static uint32_t msws() { + * x *= x; x += (w += s); return x = (x>>32) | (x<<32); + * } + * + * int main() { + * x = 0x012de1babb3c4104; + * w = 0xc8161b4202294965; + * s = 0xb5ad4eceda1ce2a9; + * for (int i=0; i<10; i++) { + * for (int j=0; j<4; j++) { + * printf("0x%08x, ", msws()); + * } + * printf("\n"); + * } + * } + */ + final long[] seed = {0x012de1babb3c4104L, 0xc8161b4202294965L, 0xb5ad4eceda1ce2a9L}; + + final int[] expectedSequence = { + 0xe7f4010b, 0x37bdb1e7, 0x05d8934f, 0x22970c75, + 0xe7432a9f, 0xd157c60f, 0x26e9b5ae, 0x3dd91250, + 0x8dbf85f1, 0x99e3aa17, 0xcb90322b, 0x29a007e2, + 0x25a431fb, 0xcc612768, 0x510db5cd, 0xeb0aec2f, + 0x05f88c18, 0xcdb79066, 0x5222c513, 0x9075045c, + 0xf11a0e0e, 0x0106ab1d, 0xe2546700, 0xdf0a7656, + 0x170e7908, 0x17a7b775, 0x98d69720, 0x74da3b78, + 0x410ea18e, 0x4f708277, 0x471853e8, 0xa2cd2587, + 0x16238d96, 0x57653154, 0x7ecbf9c8, 0xc5dd75bf, + 0x32ed82a2, 0x4700e664, 0xb0ad77c9, 0xfb87df7b, + }; + + RandomAssert.assertEquals(expectedSequence, new MiddleSquareWeylSequence(seed)); + } + + /** + * Test the self-seeding functionality. The seeding of the generator requires a high complexity + * increment for the Weyl sequence. The standard test for the statistical quality of the + * generator uses a good increment. This test ensures the generator can be created with a seed + * smaller than the seed size without exception; the statistical quality of the output is not + * assessed and expected to be poor. + */ + @Test + public void testSelfSeeding() { + // Ensure this does not throw. The output quality will be poor. + final MiddleSquareWeylSequence rng = new MiddleSquareWeylSequence(new long[0]); + rng.nextInt(); + } +}