The generator's lines
A linear congruential generator outputs the numbers xn+1 = (a·xn + c) mod m. If all consecutive pairs (xn/m, xn+1/m) are marked on a square, they fall on parallel lines. Integers s1, s2, not both zero, define the family of lines s1·u + s2·v = k + s2·c/m, k any integer, and this family covers all the pairs if s1 + a·s2 is divisible by m. The distance between neighbouring lines of the family is 1/√(s12 + s22).
The spectral test looks for the family with the sparsest lines, that is, the smallest s12 + s22: the larger it is, the denser the lines and the more evenly the pairs fill the square. Find it for each generator. The answer does not depend on c.
Input
The first line holds the number of generators q (1 ≤ q ≤ 2·104). Then q lines “a m”: 2 ≤ m ≤ 231, 1 ≤ a < m.
Output
For each generator, the smallest value of s12 + s22 over all integers s1, s2, not both zero, for which s1 + a·s2 is divisible by m.
The example holds the generators from the article: RANDU, MINSTD with the multipliers 16807 and 48271, and three generators from the scene with modulus 4096. For a = 2049 the answer is 8: the vector (−2, 2), and all the pairs lie on four lines.
Examples
Example 1
6 65539 2147483648 16807 2147483647 48271 2147483647 5 4096 2049 4096 1697 4096
2147221514 282475250 1990735345 26 8 4562
Limits
- Time · Python
- 4000 ms
- Time · Go
- 1000 ms
- Memory
- 256 MB
- Judging
- 6 tests
Walkthrough
Open it step by step: each step gives away a little more. You can jump straight to the solution, but a problem solved with one hint stays with you longer.
1Hint 1
All pairs (s1, s2) with s1 + a·s2 ≡ 0 (mod m) form a lattice in the plane. Which two vectors generate it? And how do you find the shortest nonzero vector of a lattice in a way that resembles the Euclidean algorithm?
2The idea
Any solution can be written as s2·(−a, 1) + j·(m, 0): the second term adds a multiple of m to s1. So the lattice is generated by the vectors (m, 0) and (−a, 1), and its shortest nonzero vector is needed.
In the plane it is found by Lagrange–Gauss reduction, a two-dimensional Euclidean algorithm. Keep two vectors generating the lattice: u the longer, v the shorter. Subtract from u the multiple of v nearest to the projection: u ← u − k·v, where k is (u·v)/(v·v) rounded. The lattice does not change. If the new u is shorter than v, swap their roles and continue; if not, v is the shortest. Why: after stopping |u| ≥ |v| and |u·v| ≤ v·v/2, and then any combination αu + βv with integer α ≠ 0 is no shorter than v, while with α = 0 it is a multiple of v.
3Walkthrough
Start with u = (m, 0), v = (−a, 1), swapping them if u is shorter. The rounding k = ⌊(2·u·v + v·v) / (2·v·v)⌋ is to the nearest integer. The answer is v·v at the moment of stopping.
Squared lengths and dot products are at most m2 ≤ 262 and fit in a 64-bit signed integer. But 2·u·v may no longer fit, so in Go the rounding goes through the remainder: k = ⌊d/vv⌋, plus one if the remainder is at least half of vv.
3+Solution
4Complexity
As with the Euclidean algorithm: the lengths fall geometrically, O(log m) steps per generator, O(q log m) in all. The most steps come with a close to m·(√5 − 1)/2, as with Euclid on consecutive Fibonacci numbers.
5Common mistakes
- Trying s2 up to a small bound. The shortest vector can be up to (4/3)1/4·√m long, and with m around 231 s2 reaches tens of thousands.
- Stopping when u is shorter instead of swapping. A step must always subtract the short vector from the long one, otherwise the reduction gets stuck or returns something other than the shortest.
- Rounding by dropping the fractional part. Truncation towards zero instead of to the nearest integer leaves u longer than it could be, and the loop ends too early.
Algorithms behind it
O(1) per number
2 problems
NUM
Algorithm
Euclid's algorithm
The greatest common divisor of two numbers in O(log n) divisions with remainder.
O(log min(a, b))
2 problems