The New Sudoku Players' Forum

[Red Ed]

Challenge: given a perfect(*) random number generator, can you produce a perfect random sudoku grid generator
(*: perfect = no bias, correlations etc.)

I have coded up one such generator as follows (using the Mersenne Twister as my not-quite-perfect RNG; thanks Mike Metcalf for the tip-off).

One-off work

• Define a 'stub' to be a grid with nine 1s, nine 2s and nine 3s filled in, but no other digits. There are T(3) = 5196557037312 stubs.
• Collect the stubs into isomorphism classes; two stubs are in the same class iff they are isomorphic (in the usual sense). It turns out that there are just 259272 classes, which we tabulate as follows:
  • stub(k) = a representative grid;
  • csz(k) = the class size, i.e. the number of equivalent stubs (~20 million on average);
  • ext(k) = the number of ways of extending to a full solution grid (~1.3 billion on average).
• Write the stub(.), csz(.) and ext(.) tables to a file for later reference.

Per-grid work

• Select a random class with Pr(class k) = csz(k) * ext(k) / 6670903752021072936960.
• Select a random solution number, n in {1,2,...,ext(k)}.
• Enumerate solutions to stub(k), stopping at the nth (call this solution G).
• Apply a random digit/band/stack/row/column/transpose permutation to G, and report the result.

Outline proof of unbiasedness

• Let Solns(s) be the solutions to stub s. Then for each stub s in class k, Solns(s) is isomorphic to Solns(stub(k)). All Solns(.) are distinct, so there are csz(k) * ext(k) solution grids represented by class k. (Thus, class k should be picked with probability proportional to that.)
• Since the permutation (group op) in the last step is chosen flat randomly, it can be preceded by any other group ops (chosen randomly or not) without changing the average effect. In particular, we could apply a random permutation of {1,2,3} and bands/stacks/rows/cols/transpo on stub(k) prior to enumerating solutions. This has the effect of choosing a random stub in class k without bias which, coupled with the solution-picking step, results in unbiased selection of a solution grid from the class as a whole.

The algorithm trades memory (one-off, building the isomorphism classes) for time (per grid, enumerating solutions). With the memory/time split as above, the algorithm generates a rather paltry 10 grids/sec on my machine. An obvious speed improvement would be to do more up-front work, perhaps caching some of the results on disk if the tables got too big. Classifying 4-digit templates might be one option; or trying a different stub scheme altogether (such as isomorphism classes of filled bands 1&2) might be preferable. I'd be pleased to see something better implemented.

Finally, an application. We'd like to know how many proper (solvable, not necessarily minimal) sudoku puzzles exist. To answer this, we generate grids in an unbiased manner and, for each one, try to solve a random subgrid (that is, clear each cell independently with probability .5 then ask if what's left has a unique solution). Let P be the probability that a puzzle generated this way is proper; then the number of proper puzzles in total is P * 6.67e21 * 2^81. In 178100 trials, I found 65325 proper puzzles, so a 95% confidence interval for the total number of proper puzzles is [ 5.88e45, 5.95e45 ] to 3 s.f.

[JPF]

Thanks for this interesting and opportune work.
I will try to understand what you did, if I can.

Jumping to your conclusions, it appears that my grid generation is biased related to yours : P = 0.378 against 0.367.

I have no doubt you are right.

My number generator seems not too bad.
I tested it with the estimation of Π/4 by simulating (x^2 + y^2 ≤ 1) ; [0≤x≤ 1; 0≤y≤ 1] with 10^p trials ; p=4,…,9 and more
It worked rightly.

To estimate P, I used 2 main generations :

1) The generation of patterns.
That is to pick a uniform number between 0 and 2^81-1. i.e. 81 independent binary digits (0,1).

2) The random generation of a full grid.
That is probably during this phase that I 'm faced with a bias.

JPF

[Red Ed]

Interesting point about the source of bias; yes, I guess it could be in the puzzle selection step rather than the grid generation step. I had assumed it was down to grid generation simply because the values of 'P' vary so much per solution grid, so picking those in an unbiased manner seemed important.

But your thought is definitely worth testing. Can you replace your puzzle selector with code along the following lines (just showing the relevant snippets below) ...

Code: Select all

#include "mtwist-0.6/mtwist.h"
#include "mtwist-0.6/randistrs.h"

mt_state state;
mts_goodseed(&state);

int puz[9][9];
make_random_grid(puz);   /* fill puz with a random solution grid */

unsigned long zap[3];
zap[0] = mts_lrand(&state);
zap[1] = mts_lrand(&state);
zap[2] = mts_lrand(&state);
for (i=0; i<81; i++) if (zap[i/27] & (1<<(i%27))) puz[i/9][i%9] = 0;

... and see what happens? (This example uses the code at http://www.cs.hmc.edu/~geoff/mtwist.html btw)

[Ocean]

Challenge: given a perfect(*) random number generator, can you produce a perfect random sudoku grid generator
(*: perfect = no bias, correlations etc.)

I have coded up one such generator as follows (using the Mersenne Twister as my not-quite-perfect RNG; thanks Mike Metcalf for the tip-off).

Interesting work. Maybe I will be able to understand more details of random grid generation. (Never studied the grid counting/generation algorithms in detailed detail, since this has been in good hands, and no new ideas to be expected from my side...).

meanwhile...

Number of valid puzzles* with c clues. (3)
* minimal & non minimals
[...]
Here are confidence intervals [N1, N2] (probability = 95%) for the number of valid puzzles with c clues.

JPF: I have basically replicated your numbers (interpreted as averages for a certain set of grids), with the following procedure:

1. Generate a random grid.
2. Delete randomly x clues from this grid (C=81-x).
3. If a valid puzzle, count the number of removable clues for this puzzle. (If not valid, scrap it).
4. Apply the recursive formula. (N(C-1)=r(C)*N(C)/(82-C); starting with N(81)=1, and r(C) is average number of removable clues, taken over all grids that resulted in a valid C-puzzle).

For a grid set of size about 10000: The numbers N(C) were within your 95% confidence intervals for "almost all" values of C (26<=C<=81). (Only three numbers were slightly outside the interval, one 'significant decimal' away).

I think this procedure, when lasting long enough, will produce correct average number of puzzles for the set of grids generated.

The grid generation is probably biased though: Fill in numbers "randomly" in an empty grid, until this defines a unique solution. Select this solution grid. For filling in new numbers in a multisolution grid, a random posistion is chosen, and a random valid number is chosen for that position. (Every time 'random' means: output from a pseudorandom generator.) The bias comes from the assumption that - with this procedure - some grids are generated more often than others (in the infinite case). (Can be show/proved by examples.)

Comments:
1. Interesting if the biased grid generation gives a different average number of puzzles per grid than the true average.
2. In any case, this might serve as an 'independent verification' of the procedures for finding the total number of puzzles in a certain set of grids (or the average per grid in such a set).
3. Observation: Not a single minimal puzzle was generated during this procedure.

[JPF]

Can you replace your puzzle selector with code along the following lines (just showing the relevant snippets below) ...

Unfortunately, I don’t use C ; the C-mtwist pakage can’t help me to improve my generator and I'm not in a position to make your suggested test.
Did you test this brand-new PRNG with the specific grids we looked at : SFB, SF, Pt, dukuso16, Ran11, 20-17s, RV, MC ?

I made 10^4 trials for each grid and 10^5 for dukuso16.
So it’s easy to calculate the resulting 95% confidence interval for each grid.
Is there any bias for those grids ?

JPF: I have basically replicated your numbers (interpreted as averages for a certain set of grids), with the following procedure

Well, any reason that our two programs with two different procedures give the same bias ?
My procedure is :
For a given c : N=0 ; Nv=0
1. Generate a "random grid". N=N+1
2. Delete randomly x clues from this grid (x=81-c).
3. If it's a valid puzzle, Nv=Nv+1
go to 1.
4. N(c) = Nv/N x ( 6.67 x 10 ^21) x 81 ! / (c ! x (81-c) !)

For a grid set of size about 10000: The numbers N(C) were within your 95% confidence intervals for "almost all" values of C (26<=C<=81). (Only three numbers were slightly outside the interval, one 'significant decimal' away).

What are these 3 values of c ?

JPF

[Ocean]

Well, any reason that our two programs with two different procedures give the same bias ?
My procedure is :
For a given c : N=0 ; Nv=0
1. Generate a "random grid". N=N+1
2. Delete randomly x clues from this grid (x=81-c).
3. If it's a valid puzzle, Nv=Nv+1
go to 1.
4. N(c) = Nv/N x ( 6.67 x 10 ^21) x 81 ! / (c ! x (81-c) !)

Seems that the procedures for removing clues are the same, and calculation of average should give equvalent results in the infinite case. Or more precisely: Your formula: N(c) = Nv/N x ( 6.67 x 10 ^21) x 81 ! / (c ! x (81-c) !) leads to the same result as mine: N(C-1)=r(C)*N(C)/(82-C); starting with N(81)=1 (except for the factor 6.67 x 10 ^21). The grid generation procedures probably have some bias, and this bias could be different.
Even most of the numbers are within the 95% confidence interval, the resulting P is 0.386 (versus 0.378 and 0.367).

After having 'concluded' the above, I decided to dive into details of the calculations. Result: I can not explain why it converges so 'fast'.

What are these 3 values of c ?

The point was to show convergence for a relatively small dataset. My exact numbers are not so important, but the observation that they rather quickly seem to converge to some value within your confidence intervals. The 3 c's were for c=43, c=32, c=29. Due to the nature of the recursion formula, it's expected that calculations are less accurate at lower c. It actually involves more than 40 estimated averages for r(C) multiplicated! Have no idea how to estimate accuracy in such a case. The convergence surprises me. Example: For C=27, 15 subgrids are used: {4 x r4; 1 x r5; 2 x r6; 4 x r7; 2 x r8; 1 x r9; 1 x r10 } resulting in an estimated average r(27)~6.4 reducable clues for an 'average' grid with 27 clues. For higher C, the number of subgrids increases (C28: 719/78~9.2; C29: 2309/204~11.3; C30: 5814/426~13.6; etc), and the number of reducable clues is more smoothly distributed.

Also, my maximum is at c=43, while yours is at c=42; indicating a different bias?

Anyway, here are estimated numbers r(C) and N(C). The number of decimals does not reflect the accuracy. (Did not check why machine-calculation of r(C) deviates a little bit from the hand-calculations on tagged subgrids in an output file, but they are roughly the same).

Code: Select all

The biased grid set (sample size ~44000 grids):

   C: Number of clues.
r(C): Average number of reducable clues.
N(C): Calculated average number of valid sudokus per grid.

 C   r(C)              N(C)
-----------------------------
81   81                   1
80   80                  81
79   79                3240
78   77,999593        85320
77   76,998371      1663731
76   75,995724     25620920
75   74,991198    324513397
...

35   25,290735      1,3E+22
34   23,200704      7,0E+21
33   20,972918      3,4E+21
32   18,698748      1,4E+21
31   16,410737      5,4E+20
30   13,916667      1,7E+20
29   11,606061      4,6E+19
28    9,237288      1,0E+19
27    6,538462      1,7E+18
26    -             2,1E+17


Red Ed: Would be interesting if you could provide (average/distribution of) number of reducable clues for subgrids from a grid set generated by your unbiased procedure - say for C=78 (related to unavoidable sets of size four), and C below 30 (the 'unstable' interval).

[Red Ed]

Ooh, fun! Unavoidables of size 4 seem to provide a nice way of detecting certain types of bias. There are 486 possible ways of placing a size-4 unavoidable on a grid; the average number of size-4 unavoidables is X, a number not much bigger than 10 (I am being deliberately vague). So, here's a test. Generate a zillion random grids and count the average number of size-4 unavoidables. What do you get? I'll post the exact value of X in the next few days ... anyone want to put their generator up to the challenge before then?

One of my previous attempts at a "reasonably unbiased" grid generator had an X value that was off by about 4%. gsf's program (command line: sudoku -g -f%v) is way off the mark, scoring nearly 60.

The test can easily be extended to cover different sorts of structures (not just size-4 unavoidables), which should make it quite a bit more powerful. I'll come back to that in a few days when reporting the true value of X.

( Update: to calculate the number of size-4 unavoidables in a grid, do this. Set U=0. Loop over a band B and pair of rows r1,r2 in B (9 choices); loop over a pair of stacks S1,S2 and a column c1,c2 in each (27 choices); increment U if (grid[r1,c1]==grid[r2,c2] AND grid[r1,c2]==grid[r2,c1]); increment U if (grid[c1,r1]==grid[c2,r2] AND grid[c2,r1]==grid[c1,r2]). Report the average of U over lots of grids. )

[JPF]

Unavoidables of size 4 seem to provide a nice way of detecting certain types of bias. There are 486 possible ways of placing a size-4 unavoidable on a grid; the average number of size-4 unavoidables is X, a number not much bigger than 10 (I am being deliberately vague). So, here's a test. Generate a zillion random grids and count the average number of size-4 unavoidables. What do you get? I'll post the exact value of X in the next few days ... anyone want to put their generator up to the challenge before then?

So, let's shoot first !
The Zillions are not necessary for a biased generation...

Here are my numbers with 10^6 random (biased) grids.

X is the random variable of the number of size-4 unavoidable on each grid.

E (X) = 10.325
s (X) = 4.797
Min (X) = 0
Max (X) = 36

95% confidence interval : 10.315 ≤ X ≤ 10.334
[edit : suppressed distribution of X ; to be checked ]

I'm waiting for sarcastic remarks

JPF

[coloin]

Nice idea

Here is a grid with only 3 pairs...with one common clue interestingly

Code: Select all

+---+---+---+
|137|824|596|
|695|371|428|
|248|659|713|
+---+---+---+
|351|987|642|
|476|213|985|
|982|465|137|
+---+---+---+
|563|142|879|
|719|538|264|
|824|796|351|
+---+---+---+

+---+---+---+
|...|.24|...|
|...|...|...|
|...|...|...|
+---+---+---+
|...|...|...|
|...|...|...|
|...|...|...|
+---+---+---+
|.6.|.42|8..|
|...|..8|2..|
|.2.|..6|...|
+---+---+---+

Are grids with no size-4 set really that common ?
You must mean X = 0-1 p% = 1.42

Im surprised that we havnt done this sort of analysis before - re structure of the solution grid. Perhaps it was more difficult to count up unavoidables of size greater than x=8 envolving more than 2-digits.......?

Will extraordinary grids have an extraordinary distribution of size-x unavoidables ?

C

[Red Ed]

E (X) = 10.325
s (X) = 4.797
Min (X) = 0
Max (X) = 36

Not bad -- and much better than gsf's program managed (at least with me driving it).

I'm waiting for sarcastic remarks

Ooh, wow, that's like SOOOO AMAZING with those numbers and everything. ( Will that do? )

[JPF]

Ooh, wow, that's like SOOOO AMAZING with those numbers and everything. ( Will that do? )

Perfect ! Thanks

Not bad -- and much better than gsf's program managed (at least with me driving it).

Could be better. I discovered a bug in my counting program.
I'm doing it again.

JPF

[JPF]

Hope my bug is fixed now.
So, with the same bias in the grids generation :

n = 500,000 grids
E (X) = 10.99
s (X) = 3.43

You must mean X = 0-1 p% = 1.42

It was wrong.
Now I get p(0) ~ 0.01%.

JPF

[coloin]

So this wayward result confirms a bias ?
And demonstrates the validity of Red Ed's method - await his unbiased result !

Care to see how biased dukuso's generator is ?
the program is suexg
command would be suexg 0 1000000 p>file
unfortunately this will generate puzzles [not full grids]

C

[gsf]

gsf's program (command line: sudoku -g -f%v) is way off the mark, scoring nearly 60.

if you just want a stream of grids try

Code: Select all

sudoku -g -T0x80 -f%v

[Ocean]

Illustration of a biased grid generator: Pick randomly an empty cell; Assign randomly a valid digit to that cell; (repeat until grid is filled).

"Every time (*)" this is reached:
(*) (happens less than once per life span of our galaxy on my computer)

..4.56789.57.89246689247135..3564897948723561765891423591438672432675918876912354

there are three possible resulting grids:

124356789357189246689247135213564897948723561765891423591438672432675918876912354 # p=26/168
214356789357189246689247135123564897948723561765891423591438672432675918876912354 # p=71/168
324156789157389246689247135213564897948723561765891423591438672432675918876912354 # p=71/168

The chance of getting to the first grid is ~16%, while for each of the other two ~42%.
So, it is suspected (expected) that certain grid types will be overrepresented compared to others, when using this particular method for grid generation.
#