Latin Squares in reduced form/Randomizing using Jacobson and Matthews’ Technique: Difference between revisions

From Rosetta Code
Content added Content deleted
(→‎{{header|Go}}: Added a second much quicker but less uniform version.)
(→‎{{header|Go}}: First version now far quicker thanks to advice by Nigel Galloway and second version removed as no longer necessary.)
Line 128:
 
=={{header|Go}}==
The J & M implementation is based on the C code [https://brainwagon.org/2016/05/17/code-for-generating-a-random-latin-square/ here] which has been heavily optimized following advice and clarification by Nigel Galloway (see Talk page) on the requirements of this task.
===Parts 1, 2 and (a small part of) 3===
The J & M implementation is based on the C code [https://brainwagon.org/2016/05/17/code-for-generating-a-random-latin-square/ here]. Although this produces satisfactory results for Parts 1 and 2, it is glacial for Part 3 and could only generate 5 order 42 Latin squares in around 22 seconds. Optimizing the code somewhat gets this down to about 14 seconds which suggests that the full 750 could be generated in about 35 minutes though Part 4 is still well out of reach.
 
AllPart timings4 areis fortaking about 6.5 seconds on my Celeron @1.6 GHz andbut will therefore be much faster on a more modern machine. thoughBeing probablyable stillto farcompute toorandom, slowuniformly fordistributed, Latin squares of order 256 reasonably quickly is interesting from a secure communications or cryptographic standpoint as the latersymbols Partsof such a square can represent the 256 characters of the various extended ASCII encodings.
<lang go>package main
 
Line 233 ⟶ 232:
n := len(c[0])
proper := true
itersvar :=rx, nry, *rz n * nint
for i := 0; i < iters; i++ {
var rx, ry, rz= intrand.Intn(n)
forry {= rand.Intn(n)
rxrz = rand.Intn(n)
if c[rx][ry][rz] == 0 rand.Intn(n){
rz = rand.Intn(n)break
if c[rx][ry][rz] == 0 {}
}
for {
var ox, oy, oz int
for ; ox < n; ox++ {
if c[ox][ry][rz] == 1 {
break
}
}
forif !proper && rand.Intn(2) == 0 {
varfor ox,++; oy,ox oz< n; ox++ int{
for ; ox < n; ox++ {
if c[ox][ry][rz] == 1 {
break
}
}
}
if !proper && rand.Intn(2) == 0 {
for ox++; ox < n; ox++ {
if c[ox][ry][rz] == 1 {
break
}
}
}
 
for ; oy < n; oy++ {
if c[rx][oy][rz] == 1 {
break
}
}
if !proper && rand.Intn(2) == 0 {
for oy++; oy < n; oy++ {
if c[rx][oy][rz] == 1 {
break
}
}
}
if !proper && rand.Intn(2) == 0 {
 
for oy++; oy < n; oy++ {
for ; oz < n; if c[rx][oy][rz] == 1oz++ {
if c[rx][ry][oz] == 1 break{
}break
}
}
}
 
if !proper && rand.Intn(2) for== ; oz < n; oz++0 {
for oz++; oz < n; oz++ {
if c[rx][ry][oz] == 1 {
break
}
}
}
if !proper && rand.Intn(2) == 0 {
for oz++; oz < n; oz++ {
if c[rx][ry][oz] == 1 {
break
}
}
}
 
c[rx][ry][rz]++
c[rx][oy][oz]++
c[ox][ry][oz]++
c[ox][oy][rz]++
 
c[rx][ry][oz]--
c[rx][oy][rz]--
c[ox][ry][rz]--
c[ox][oy][oz]--
 
if c[ox][oy][oz] < 0 {
rx, ry, rz = ox, oy, oz
proper = false
} else {
proper = true
break
}
}
}
Line 391 ⟶ 387:
 
// part 3
fmt.Println("\nPART 3: 5750 latin squares of order 42, showing the last one:\n")
start := time.Now()
var m42 matrix
c = makeCube(nil, 42)
for i := 1; i <= 5750; i++ {
shuffleCube(c)
m42if i == 750 toMatrix(c){
m42 = toMatrix(c)
}
}
elapsed := time.Since(start)
printMatrix(m42)
 
fmt.Printf("Took %s\n", elapsed)
// part 4
fmt.Println("\nPART 4: 1000 latin squares of order 256:\n")
start := time.Now()
c = makeCube(nil, 256)
for i := 1; i <= 1000; i++ {
shuffleCube(c)
}
elapsed := time.Since(start)
fmt.Printf("Generated in %s\n", elapsed)
}</lang>
 
Line 408 ⟶ 413:
<pre>
PART 1: 10,000 latin Squares of order 4 in reduced form:
 
1 2 3 4
2 1 4 3
3 4 1 2
4 3 2 1
 
Occurs 2509 times
 
1 2 3 4
2 4 1 3
3 1 4 2
4 3 2 1
 
Occurs 2571 times
 
1 2 3 4
2 3 4 1
3 4 1 2
4 1 2 3
 
Occurs 2446 times
 
1 2 3 4
Line 435 ⟶ 419:
4 3 1 2
 
Occurs 24742550 times
 
 
PART 2: 10,000 latin squares of order 5 in reduced form:
 
1(183), 2(189), 3(162), 4(188), 5(191), 6(171), 7(154), 8(212),
9(169), 10(177), 11(198), 12(174), 13(181), 14(177), 15(189), 16(162),
17(202), 18(188), 19(167), 20(194), 21(167), 22(158), 23(157), 24(154),
25(183), 26(160), 27(192), 28(186), 29(187), 30(171), 31(159), 32(182),
33(156), 34(198), 35(165), 36(170), 37(166), 38(163), 39(182), 40(161),
41(198), 42(188), 43(176), 44(190), 45(217), 46(178), 47(153), 48(199),
49(178), 50(158), 51(192), 52(205), 53(177), 54(184), 55(183), 56(179)
 
 
PART 3: 5 latin squares of order 42, showing the last one:
 
18 21 11 26 2 8 25 39 38 24 3 19 33 37 5 16 30 4 9 40 32 31 36 22 1 42 7 12 27 35 23 13 17 14 15 10 41 29 20 34 6 28
41 30 9 40 37 11 16 14 2 32 42 33 7 5 12 1 18 34 23 19 26 22 35 25 36 24 4 21 28 6 39 15 20 8 27 38 17 3 13 10 31 29
8 29 34 2 19 3 6 4 23 30 36 38 11 39 26 28 10 7 22 24 42 25 16 32 17 14 5 9 33 41 31 40 18 27 13 1 37 21 35 12 20 15
30 15 8 18 28 7 31 2 16 10 38 4 26 29 17 25 9 35 14 20 41 33 39 11 24 5 3 36 21 27 1 23 37 40 42 22 12 19 34 32 13 6
39 41 20 17 27 19 24 38 34 16 6 32 29 13 28 3 25 8 1 11 31 12 10 33 23 21 37 2 40 5 18 22 14 7 26 35 9 42 15 36 4 30
2 1 23 15 26 35 36 29 13 6 27 31 12 30 40 7 39 19 3 34 21 17 42 20 14 38 9 8 41 11 22 33 4 5 24 28 18 37 25 16 32 10
28 31 7 10 40 24 30 11 35 23 33 17 32 18 22 12 42 25 13 15 8 21 9 16 41 3 27 5 38 36 6 20 19 4 1 29 39 34 26 2 14 37
6 26 3 12 13 40 38 30 33 17 20 29 28 9 34 15 24 32 39 14 11 8 25 21 19 1 35 23 5 16 2 41 36 31 7 37 22 4 10 18 42 27
19 5 26 31 36 12 33 6 4 1 21 27 38 17 7 20 34 14 32 22 3 2 8 10 9 35 16 18 39 15 37 42 41 29 28 24 30 13 11 40 23 25
5 22 36 11 29 42 41 17 39 40 4 9 8 25 1 21 3 24 7 23 35 13 27 26 38 12 33 31 14 28 15 18 2 34 32 20 10 30 6 19 37 16
33 36 32 28 11 17 14 9 25 21 23 7 15 3 24 2 4 22 31 26 30 10 34 41 27 19 38 29 1 40 20 35 13 18 6 39 8 16 42 37 12 5
16 25 10 3 23 26 20 32 17 5 31 40 36 1 4 27 14 39 11 28 12 9 24 19 33 8 42 22 2 30 38 21 34 15 29 7 35 6 37 13 41 18
40 28 17 16 39 25 37 41 22 15 30 8 42 20 11 19 35 26 10 29 7 5 13 18 21 6 36 32 12 24 27 4 33 23 9 3 38 2 31 14 1 34
25 37 1 29 7 28 9 19 26 22 13 41 21 2 6 39 12 3 16 27 17 42 14 36 18 15 30 4 11 20 35 8 24 32 10 34 40 31 33 38 5 23
37 34 18 33 35 38 27 23 41 8 28 10 39 40 31 5 11 1 30 9 24 29 3 7 25 32 6 16 26 13 14 19 42 17 2 15 20 36 12 4 22 21
15 3 37 35 38 21 32 42 36 34 16 18 30 33 19 41 20 5 27 25 13 7 28 1 2 11 24 14 6 17 8 10 40 22 39 12 4 9 23 29 26 31
10 35 13 1 42 6 34 20 9 4 5 14 22 23 38 30 7 11 28 21 27 32 2 15 8 39 29 40 37 26 25 24 12 19 41 18 3 17 16 31 36 33
23 33 39 20 25 9 40 34 30 35 15 42 37 14 3 18 21 2 8 12 36 4 17 5 26 31 10 6 13 1 24 16 22 41 38 11 27 28 19 7 29 32
12 27 16 39 30 4 18 15 14 20 35 1 40 8 9 23 38 36 41 2 34 11 19 6 29 17 22 10 25 37 13 28 31 33 5 42 26 7 32 21 24 3
24 9 38 42 41 34 21 36 18 25 8 16 1 12 10 17 27 37 19 32 20 26 33 13 30 4 28 15 29 39 7 31 5 3 35 23 6 22 14 11 40 2
22 39 19 13 9 23 11 7 20 28 10 5 18 36 41 14 40 27 15 16 38 6 31 34 3 2 17 37 35 33 30 12 29 42 25 4 24 32 21 1 8 26
9 40 27 36 10 20 4 1 42 14 2 13 5 41 37 8 15 18 25 7 19 23 6 35 32 29 26 3 22 12 33 30 28 21 16 31 34 11 17 24 38 39
3 42 25 24 8 5 17 22 1 13 14 36 20 35 2 31 26 29 37 33 18 34 11 30 15 40 19 27 7 10 28 32 9 16 4 6 23 12 38 39 21 41
20 16 29 9 1 2 15 5 31 18 19 37 35 26 32 24 6 30 4 42 22 28 23 3 12 13 25 39 36 34 41 14 10 38 40 21 11 27 7 8 33 17
17 32 28 41 6 30 8 37 21 29 24 34 10 38 18 35 19 16 33 31 25 14 15 40 4 9 12 42 23 7 3 5 27 26 22 2 36 1 39 20 11 13
11 4 5 22 16 15 12 25 40 33 1 39 31 7 23 42 36 10 29 6 14 27 30 38 28 37 13 26 34 9 32 3 35 20 21 41 2 24 18 17 19 8
32 7 35 6 22 13 28 16 37 26 17 3 24 15 25 40 33 9 18 5 10 20 12 4 11 41 14 30 8 21 29 34 38 2 31 36 1 23 27 42 39 19
34 18 22 21 32 14 39 26 7 38 12 24 27 10 36 13 29 15 17 30 33 1 41 37 42 25 20 11 19 23 4 2 16 6 8 5 31 40 3 28 35 9
31 23 40 27 3 41 29 24 11 39 25 21 34 32 15 38 13 6 12 10 9 16 37 8 20 28 1 17 4 18 5 26 7 36 33 19 42 35 2 22 30 14
38 8 33 19 15 10 23 12 32 42 7 2 6 4 14 36 17 31 35 18 37 39 1 29 16 22 21 13 24 25 34 27 3 11 20 9 5 41 30 26 28 40
13 20 42 4 17 32 3 28 15 2 37 22 14 16 29 11 1 21 26 39 6 24 5 31 40 27 34 41 9 38 12 7 30 25 19 8 33 10 36 23 18 35
27 2 30 37 31 22 35 8 19 36 41 12 4 34 39 29 16 23 38 13 28 40 21 24 6 26 15 7 20 42 9 17 25 10 18 14 32 5 1 33 3 11
1 17 24 32 20 36 10 31 29 37 22 28 23 11 27 34 41 33 42 38 4 3 40 2 5 7 39 35 15 19 26 9 6 13 14 30 21 18 8 25 16 12
26 14 15 30 24 18 42 10 8 7 32 25 3 19 33 6 31 28 40 36 5 35 4 23 34 20 2 38 16 29 21 1 11 12 37 27 13 39 41 9 17 22
35 38 41 23 4 39 26 13 28 3 18 20 19 6 16 32 2 12 5 8 1 37 22 27 10 34 40 33 17 31 42 29 21 30 11 25 7 14 24 15 9 36
29 24 6 7 21 16 22 33 12 31 11 23 13 27 8 26 28 38 34 35 40 41 18 14 39 10 32 1 3 4 19 37 15 9 36 17 25 20 5 30 2 42
36 12 31 14 18 29 5 40 10 27 39 30 17 24 21 22 37 13 6 3 15 38 32 9 7 23 8 28 42 2 16 11 1 35 34 33 19 26 4 41 25 20
14 11 4 25 34 1 7 27 5 19 9 26 41 42 30 33 22 17 21 37 2 15 20 12 35 16 18 24 31 3 40 36 32 39 23 13 29 8 28 6 10 38
4 19 14 5 33 37 2 35 24 41 29 15 9 21 20 10 23 42 36 17 39 18 38 28 13 30 31 34 32 8 11 6 26 1 12 40 16 25 22 3 27 7
7 6 21 38 14 33 13 18 3 12 40 11 2 22 35 37 5 41 20 4 29 19 26 42 31 36 23 25 10 32 17 39 8 24 30 16 28 15 9 27 34 1
21 13 2 8 12 27 19 3 6 11 34 35 25 31 42 9 32 40 24 1 16 36 7 39 37 18 41 20 30 22 10 38 23 28 17 26 14 33 29 5 15 4
42 10 12 34 5 31 1 21 27 9 26 6 16 28 13 4 8 20 2 41 23 30 29 17 22 33 11 19 18 14 36 25 39 37 3 32 15 38 40 35 7 24
 
Took 13.714977269s
</pre>
 
===Parts 1, 2, 3 and (a small part of) 4===
The J & M implementation this time is based on the optimized Go code [https://blog.paulhankin.net/latinsquares/ here] which is more than 100 times quicker than before. The full 750 order 42 latin squares are now generated in around 17 seconds. Even Part 4 is now within reach with 10 order 256 latin squares being produced in 70 seconds suggesting the full 1,000 could be produced in under 2 hours.
 
However, this is at the expense of Part 1 now producing non-uniform results and Part 2 (though more uniform) being less than ideal. This doesn't matter so much for Parts 3 and 4 where the number of possibilities is astronomic in any case though it is still, of course, a concern.
<lang go>package main
 
import (
"fmt"
"math/rand"
"time"
)
 
type (
vector []int
matrix []vector
)
 
func rand2(a, b int) (int, int) {
if rand.Intn(2) == 0 {
return a, b
}
return b, a
}
 
func toReduced(m matrix) matrix {
n := len(m)
r := make(matrix, n)
for i := 0; i < n; i++ {
r[i] = make(vector, n)
copy(r[i], m[i])
}
for j := 0; j < n-1; j++ {
if r[0][j] != j {
for k := j + 1; k < n; k++ {
if r[0][k] == j {
for i := 0; i < n; i++ {
r[i][j], r[i][k] = r[i][k], r[i][j]
}
break
}
}
}
}
for i := 1; i < n-1; i++ {
if r[i][0] != i {
for k := i + 1; k < n; k++ {
if r[k][0] == i {
for j := 0; j < n; j++ {
r[i][j], r[k][j] = r[k][j], r[i][j]
}
break
}
}
}
}
return r
}
 
// 'm' is assumed to be 0 based
func printMatrix(m matrix) {
n := len(m)
for i := 0; i < n; i++ {
for j := 0; j < n; j++ {
fmt.Printf("%2d ", m[i][j]+1) // back to 1 based
}
fmt.Println()
}
fmt.Println()
}
 
// converts 4 x 4 matrix to 'flat' array
func asArray16(m matrix) [16]int {
var a [16]int
k := 0
for i := 0; i < 4; i++ {
for j := 0; j < 4; j++ {
a[k] = m[i][j]
k++
}
}
return a
}
 
// converts 5 x 5 matrix to 'flat' array
func asArray25(m matrix) [25]int {
var a [25]int
k := 0
for i := 0; i < 5; i++ {
for j := 0; j < 5; j++ {
a[k] = m[i][j]
k++
}
}
return a
}
 
// 'a' is assumed to be 0 based
func printArray16(a [16]int) {
for i := 0; i < 4; i++ {
for j := 0; j < 4; j++ {
k := i*4 + j
fmt.Printf("%2d ", a[k]+1) // back to 1 based
}
fmt.Println()
}
fmt.Println()
}
 
func shuffleMatrix(xy matrix) {
n := len(xy)
xz := make(matrix, n)
yz := make(matrix, n)
for i := 0; i < n; i++ {
yz[i] = make([]int, n)
xz[i] = make([]int, n)
}
for i := 0; i < n; i++ {
for j := 0; j < n; j++ {
k := xy[i][j]
xz[i][k] = j
yz[j][k] = i
}
}
 
var mxy, mxz, myz int
var m [3]int
proper := true
minIter := n * n * n
for iter := 0; iter < minIter || !proper; iter++ {
var i, j, k, i2, j2, k2 int
var i2_, j2_, k2_ int
 
if proper {
// Pick a random 0 in the array
i, j, k = rand.Intn(n), rand.Intn(n), rand.Intn(n)
for xy[i][j] == k {
i, j, k = rand.Intn(n), rand.Intn(n), rand.Intn(n)
}
// find i2 such that [i2, j, k] is 1. same for j2, k2
i2 = yz[j][k]
j2 = xz[i][k]
k2 = xy[i][j]
i2_, j2_, k2_ = i, j, k
} else {
i, j, k = m[0], m[1], m[2]
// find one such that [i2, j, k] is 1, same for j2, k2.
// each is either the value stored in the corresponding
// slice, or one of our three temporary "extra" 1s.
// That's because (i, j, k) is -1.
i2, i2_ = rand2(yz[j][k], myz)
j2, j2_ = rand2(xz[i][k], mxz)
k2, k2_ = rand2(xy[i][j], mxy)
}
 
proper = xy[i2][j2] == k2
if !proper {
m = [3]int{i2, j2, k2}
mxy = xy[i2][j2]
myz = yz[j2][k2]
mxz = xz[i2][k2]
}
 
xy[i][j] = k2_
xy[i][j2] = k2
xy[i2][j] = k2
xy[i2][j2] = k
 
yz[j][k] = i2_
yz[j][k2] = i2
yz[j2][k] = i2
yz[j2][k2] = i
 
xz[i][k] = j2_
xz[i][k2] = j2
xz[i2][k] = j2
xz[i2][k2] = j
}
}
 
// 'from' matrix is assumed to be 1 based
func makeMatrix(from matrix, n int) matrix {
m := make(matrix, n)
for i := 0; i < n; i++ {
m[i] = make(vector, n)
for j := 0; j < n; j++ {
var k int
if from == nil {
k = (i + j) % n
} else {
k = from[i][j] - 1
}
m[i][j] = k
}
}
return m
}
 
func main() {
rand.Seed(time.Now().UnixNano())
 
// part 1
fmt.Println("PART 1: 10,000 latin Squares of order 4 in reduced form:\n")
from := matrix{{1, 2, 3, 4}, {2, 1, 4, 3}, {3, 4, 1, 2}, {4, 3, 2, 1}}
freqs4 := make(map[[16]int]int, 10000)
m := makeMatrix(from, 4)
for i := 1; i <= 10000; i++ {
shuffleMatrix(m)
rm := toReduced(m)
a16 := asArray16(rm)
freqs4[a16]++
}
for a, freq := range freqs4 {
printArray16(a)
fmt.Printf("Occurs %d times\n\n", freq)
}
 
// part 2
fmt.Println("\nPART 2: 10,000 latin squares of order 5 in reduced form:")
from = matrix{{1, 2, 3, 4, 5}, {2, 3, 4, 5, 1}, {3, 4, 5, 1, 2},
{4, 5, 1, 2, 3}, {5, 1, 2, 3, 4}}
freqs5 := make(map[[25]int]int, 10000)
m = makeMatrix(from, 5)
for i := 1; i <= 10000; i++ {
shuffleMatrix(m)
rm := toReduced(m)
a25 := asArray25(rm)
freqs5[a25]++
}
count := 0
for _, freq := range freqs5 {
count++
if count > 1 {
fmt.Print(", ")
}
if (count-1)%8 == 0 {
fmt.Println()
}
fmt.Printf("%2d(%3d)", count, freq)
}
fmt.Println("\n")
 
// part 3
fmt.Println("\nPART 3: 750 latin squares of order 42, showing the last one:\n")
start := time.Now()
m = makeMatrix(nil, 42)
for i := 1; i <= 750; i++ {
shuffleMatrix(m)
}
elapsed := time.Since(start)
printMatrix(m)
fmt.Printf("Took %s\n", elapsed)
 
// part 4
fmt.Println("\nPART 5: 10 latin squares of order 256:\n")
start = time.Now()
m = makeMatrix(nil, 256)
for i := 1; i <= 10; i++ {
shuffleMatrix(m)
}
elapsed = time.Since(start)
fmt.Printf("Took %s\n", elapsed)
}</lang>
 
{{out}}
Sample run:
<pre>
PART 1: 10,000 latin Squares of order 4 in reduced form:
 
1 2 3 4
Line 770 ⟶ 426:
4 3 2 1
 
Occurs 29762430 times
 
1 2 3 4
2 1 4 3
3 4 21 12
4 3 12 21
 
Occurs 30902494 times
 
1 2 3 4
Line 784 ⟶ 440:
4 1 2 3
 
Occurs 30532526 times
 
1 2 3 4
2 1 4 3
3 4 1 2
4 3 2 1
 
Occurs 881 times
 
 
PART 2: 10,000 latin squares of order 5 in reduced form:
 
1(237165), 2(170173), 3(160167), 4(181204), 5(174173), 6(211165), 7(191215), 8(153218),
9(162168), 10(153157), 11(175205), 12(184152), 13(169187), 14(162173), 15(174215), 16(168185),
17(177179), 18(178176), 19(185179), 20(169160), 21(234150), 22(186166), 23(173191), 24(180181),
25(184179), 26(166192), 27(176187), 28(172186), 29(154176), 30(145196), 31(148141), 32(165187),
33(166165), 34(175189), 35(154147), 36(176175), 37(166172), 38(162), 39(180), 40(206172),
41(191189), 42(236159), 43(179197), 44(223158), 45(165178), 46(148179), 47(181193), 48(173175),
49(164207), 50(226174), 51(175181), 52(177179), 53(239193), 54(169171), 55(171153), 56(182204)
 
 
PART 3: 750 latin squares of order 42, showing the last one:
 
829 23 132 39 4217 41 34 330 15 268 33 939 18 117 3020 1927 12 7 106 2231 14 540 2935 3425 21 129 1710 2032 3619 16 224 2842 43 4026 33 315 3523 37 141 28 1 274 13 38 2518 2421 37 622 3215 1636 11
3617 1615 2311 3031 10 349 3138 26 610 51 1128 3937 27 148 3734 2641 2521 22 12 38 175 35 136 413 20 3529 42 218 3 18 2119 24 1539 32 727 1323 3216 4025 33 8 194 4240 96 28 412 3330 22 297 14
2336 1042 35 139 37 24 1815 34 3537 1618 2132 1125 3322 31 5 12 404 17 36 273 19 613 4211 25 398 2023 12 924 2928 1527 1416 71 86 26 9 329 4140 22 327 25 19 312 1314 30 26 441 2810 21 33 38 20
421 13 916 3542 20 253 1032 13 322 26 27 39 1417 15 1211 4225 34 1937 29 6 2319 3010 1812 36 7 131 18 336 79 3839 2441 2830 2140 35 833 1622 40 331 4128 3138 2624 58 34 223 22 114 3720 1714 5
122 39 13 40 147 2638 59 1134 2141 28 29 2737 36 35 4 256 21 626 2217 1516 18 4 30 3540 20 9 108 3915 3325 19 1632 32 2011 28 723 3424 31 210 1242 23 313 3727 1712 33 814 41 321 2429 38 425 18
33 3836 34 21 153 4213 40 194 29 7 414 2 3529 24 286 1012 2731 37 1723 26 17 9 8 320 3132 1321 11 1419 41 2337 30 395 3638 30 525 11 624 2535 42 127 18 16 2039 15 710 22 828 18 321 12 229 40
13 14 2931 16 377 3022 2539 4123 3932 3134 3816 33 524 4 40 42 212 2225 3435 2326 3218 28 2711 73 2415 21 20 9 3513 1819 1 10 8 2 341 3329 12 10 26 196 17 30 4 5 638 2137 15 428 2027 36 11
59 41 173 24 286 30 719 39 14 3616 14 1315 29 828 1623 24 932 10 18 42 26 20 341 37 38 240 3234 29 118 3425 38 192 22 40 31 65 2517 2726 36 433 3013 21 12 1535 10 237 3920 3511 3327 2142 1
11 2 2218 28 27 355 15 196 87 4140 2435 33 313 20 26 368 2334 742 39 3437 1833 1326 23 122 3813 2514 54 12 4 3215 17 25 336 4031 4216 3029 38 919 3732 1041 14 161 2127 2924 11 230 69 1210 21
3827 2034 2419 15 233 2322 35 2736 10 189 4130 3214 17 221 24 7 318 3038 2842 3641 3539 19 7 6 4240 4 1337 4011 3423 3929 2526 18 12 21 3 921 1435 16 120 1510 31 525 1617 3728 26 336 2932 11 2 813
4241 3716 15 311 35 722 13 20 1829 96 38 4 285 1024 19 310 1625 27 617 1318 11 832 59 33 407 12 242 1136 27 414 34 140 2321 33 212 38 328 1430 2115 1742 3537 3923 2514 2926 30 193 2639 2231 34 3628
41 7 18 221 1715 3816 4027 2931 3018 1024 3320 78 2536 2138 10 34 9 12 114 42 429 32 24 143 2726 39 2 5 122 3541 4221 2637 1630 1914 2311 3933 35 825 3123 1340 28 3213 3619 2017 37 156 32 5 612
17 1 3310 1820 1032 2123 11 5 430 2712 68 23 419 2821 3736 15 1914 1318 4037 33 731 1626 39 2941 3216 30 126 24 822 2635 129 42 3827 28 2 223 3638 3411 92 35 317 34 3 4 540 2019 1417 2513 2425
28 6 1732 42 211 20 8 3340 27 25 5 3441 22 1517 2316 1326 29 115 35 417 2923 4236 1939 1234 28 613 1618 38 3 10 3937 20 8 14 40 11 362 31 26 4 2124 18 325 25 2419 9 21 338 71 33 12 30 3735
1235 3940 3630 35 819 21 2812 1317 2 64 22 23 27 3 520 11 31 299 15 8 23 24 42 14 1710 3339 3228 3026 4229 2433 3413 41 916 2634 25 1932 37 41 7 418 40 105 16 3815 72 1636 2038 18 1 31
1615 25 39 326 40 61 28 820 29 2021 17 12 105 3513 18 30 1722 4110 13 378 32 263 1925 15 346 72 3117 2236 38 431 3314 2919 18 2435 23 4212 27 11 39 524 28 144 3841 3632 2129 34 942 16 37 33
25 3 40 206 3826 2912 1432 24 261 1113 42 358 42 37 39 13 225 7 39 2316 31 1035 5 1629 21 1724 27 634 17 914 12 27 412 15 18 1911 28 33 20 38 18 22 39 140 23 410 34 3631 30 41 836 3219 4
2431 1538 2536 3321 1116 26 128 3530 15 3 3732 41 818 20 341 26 146 2829 29 3017 22 425 2935 31 197 4140 3227 2337 4013 3820 23 622 1611 1019 1712 42 934 39 218 1810 14 425 1239 3624 13 4 533 7 272
1440 29 4 522 4238 3135 3711 21 2317 31 9 271 1728 3919 37 7 322 42 824 3514 4112 13 130 1133 25 2234 3032 27 36 1639 23 33 9 615 1310 18 19 208 15 245 40 266 3841 3426 1216 29 4 107 20 3 2823
40 5 2617 3139 32 4 126 1914 31 237 35 711 3538 3 91 1830 19 636 2320 33 15 28 38 16 21 3729 30 149 86 2925 27 4 362 4213 41 534 2724 2212 1110 2032 1322 12 337 2428 4118 1040 3442 1723 39 258
21 2 8 29 2024 26 931 2221 1839 3823 11 14 1619 10 120 3315 27 377 35 532 2638 12 241 12 325 4022 2816 19 154 32 256 1740 3142 41 618 4130 28 4 132 1117 36 39 343 2313 3037 3533 4227 10 5 734 9
1011 2125 14 617 3618 24 219 3932 3733 4231 87 2026 15 242 3821 4120 30 115 3227 23 441 2529 2735 4039 3328 2334 1612 3110 94 11 198 1842 34 5 13 37 7 9 516 1740 22 301 1236 1438 35 283 26 296 322
726 1121 3818 2825 1929 3215 36 141 3113 19 5 2 34 2623 38 327 2241 33 403 1710 22 217 84 16 1211 42 412 24 108 21 136 18 295 2335 30 2739 2537 2014 69 4124 3736 1533 4220 97 31 128 3540 3932
3025 1927 12 433 3417 1235 2224 10 339 28 710 42 221 38 13 9 312 15 534 16 20 3 818 11 5 1 2831 41 25 7 623 37 134 35 361 32 426 2922 3814 1819 2736 2440 15 3937 26 1738 2130 4032 1420 2311 39 29
223 3519 1225 19 2730 2637 3238 2240 1314 41 731 2517 47 29 174 3916 11 9 181 31 366 2133 15 5 824 23 38 102 3 30 378 21 5 2029 34 32 28 22 615 4220 2412 1135 3318 1636 1439 1927 4010 4113 26 42
1934 31 169 10 713 41 282 38 116 2422 4031 4226 1240 36 1 814 2741 37 393 3511 12 437 1332 2627 29 935 19 5 2330 33 28 638 2021 25 1 7 2 325 3016 21 298 3436 1415 20 342 23 17 2239 1518 25 184 1024
9 1220 11 1837 1328 2541 30 288 10 315 3236 4012 20 2326 33 2439 3632 1913 26 381 25 5 399 4142 19 7 3 8 226 37 1524 14 15 4223 67 2927 1638 10 342 2130 4 22 34 35 31 18 229 2716 1740 21 417
26 28 1030 21 923 24 429 33 3 6 171 2310 34 6 533 38 252 3627 2940 3114 3534 2131 2215 4119 4237 18 2 159 40 304 2013 1635 39 248 12 320 1236 1116 1417 1932 1341 7 1825 3239 42 8 375 2726 22 111 38
31 32 3712 11 308 2940 11 16 2023 2128 2518 2442 4041 1730 38 353 38 433 92 1522 3319 36 104 25 637 12 281 2631 2720 1336 19 225 1 149 7 13 817 14 5 6 3 4127 39 1834 2324 3435 4221 15 226 29 10
318 3637 41 910 1536 1728 11 442 3913 3134 40 222 35 1 5 7 22 240 2139 38 3 630 16 131 38 5 3327 20 3417 2819 1433 4126 3215 3525 23 106 2421 29 823 4231 12 294 27 309 2532 37 118 1812 1914 2624 16
1839 24 829 2737 2225 1619 3833 27 117 2916 3210 2640 3036 4212 2830 41 11 14 21 64 34 13 15 35 2 5 3332 19 241 3931 3714 4138 1718 42 4 233 10 9 7 20 6 920 3621 4022 38 2513 1223 35 228 3126
3419 2414 30 5 4 8 540 31 3 729 12 366 3821 1926 23 615 4216 33 128 3231 3338 1413 41 109 17 1127 2612 3710 2211 16 297 2524 20 8 201 39 354 41 239 1825 2330 2122 2832 13 402 2736 42 335 34 918 1537
37 7 4132 34 5 188 1236 2041 25 142 3512 2124 3216 1339 1633 31 4 1513 10 246 3928 11 2738 22 20 342 3640 2818 31 339 3410 14 929 1726 29 381 3023 2615 4221 4027 19 217 11 1 3 6 5 825 2335 1930
27 4 41 427 26 412 14 2442 17 3715 3438 30 3635 12 825 3213 4028 1239 4220 25 10 181 3816 2833 1336 3123 67 1540 37 732 2924 2110 2531 98 33 6 121 14 326 29 11 20 223 19 399 18 519 2322 16 3534
3938 2735 42 2523 36 84 10 912 11 4 155 1821 1327 1432 3417 2025 3024 16 2118 28 40 720 22 6 242 4014 2922 30 26 3739 1733 41 128 1137 57 3213 3134 10 1 329 2315 3519 62 2441 19 389 31 116 33 3
2930 33 531 3324 12 41 336 3519 23 4032 42 194 37 629 11 2 1534 39 2516 1814 3221 3042 41 6 426 37 211 1438 24 273 2817 1122 10 132 1640 18 120 22 369 3835 1728 13 8 205 27 715 3125 10 9 268 34 7
2242 3028 32 403 14 6 231 3325 16 1222 1734 3723 2139 19 299 1835 35 3440 26 9 1436 7 3610 3531 2632 4221 2513 1041 30 5 1518 4 38 41 246 37 229 2817 33 812 1311 2720 1119 24 1 398 31 202 15 27
2016 42 5 338 19 6 910 1727 12 384 30 373 2840 2918 4111 2413 1122 3935 1 821 2 3234 2336 78 4023 2130 3417 39 542 27 267 1415 3537 3632 1620 1526 2531 19 633 1828 2229 33 109 3125 14 424 1312 41
3524 23 333 18 714 62 3225 1639 1529 19 5 199 12 295 4128 20 827 3138 20 147 33 428 2531 2311 15 110 1835 2234 2712 3916 2132 1017 21 36 3040 28 37 343 26 30 242 38 1 9 1122 4 17 13 2437 4041 6
1512 20 1 192 2329 34 365 33 42 7 24 33 104 2618 22 1114 19 936 21 389 2527 4037 2028 3126 30 838 23 310 1811 431 17 3034 15 713 3541 3740 27 283 39 1 5 8 16 29 6 235 32 21 625 12 41 13 1439
3213 68 21 139 3927 37 242 26 6 120 1725 1639 3140 3029 1032 18 4 5 328 30 24 2741 14 1922 2033 3421 3735 25 181 1215 16 836 10 9 224 3534 11 15 3326 38 11 7 362 23 423 2912 4031 28 417 23 517 19
6 34 14 2610 22 13 414 3920 25 367 18 435 19 185 1038 2313 30 142 24 206 36 923 32 821 38 172 4229 24 5 111 12 318 33 226 40 28 3727 4039 19 3 329 3525 11 731 3317 2914 2737 16 2112 1534 3041 3 15
 
Took 16.74757144s
 
PART 54: 101000 latin squares of order 256:
 
Generated in 6.581088256s
Took 1m9.236940447s
</pre>

Revision as of 09:47, 15 August 2019

Latin Squares in reduced form/Randomizing using Jacobson and Matthews’ Technique is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

Section 3.3 of [Generalised 2-designs with Block Size 3(Andy L. Drizen)] describes a method of generating Latin Squares of order n attributed to Jacobson and Matthews. The purpose of this task is to produce a function which given a valid Latin Square transforms it to another using this method.

part 1

Use one of the 4 Latin Squares in reduced form of order 4 as X0 to generate 10000 Latin Squares using X(n-1) to generate X(n). Convert the resulting Latin Squares to their reduced form, display them and the number of times each is produced.

part 2

As above for order 5, but do not display the squares. Generate the 56 Latin Squares in reduced form of order 5, confirm that all 56 are produced by the Jacobson and Matthews technique and display the number of each produced.

part 3

Generate 750 Latin Squares of order 42 and display the 750th.

part 4

Generate 1000 Latin Squares of order 256. Don't display anything but confirm the approximate time taken and anything else you may find interesting

F#

The Functions

<lang fsharp> // Jacobson and Matthews technique for generating Latin Squares. Nigel Galloway: August 5th., 2019 let R=let N=System.Random() in (fun n->N.Next(n)) let jmLS α X0= 

 let X0=Array2D.copy X0
 let N=let N=[|[0..α-1];[α-1..(-1)..0]|] in (fun()->N.[R 2])
 let rec randLS i j z i' j' z'=
   X0.[i,j']<-z'; X0.[i',j]<-z'
   match X0.[i',j']=z' with
    true->X0.[i',j']<-z; X0
   |false->randLS i' j' z' (List.find(fun n->X0.[n,j']=z')(N())) (List.find(fun n->X0.[i',n]=z')(N())) (if (R 2)=0 then let t=X0.[i',j'] in X0.[i',j']<-z; t else z)
 let i,j=R α,R α
 let z  =let z=1+(R (α-1)) in if z<X0.[i,j] then z else 1+(z+1)%α
 let i',j',z'=let N=[0..α-1] in (List.find(fun n->X0.[n,j]=z) N,List.find(fun n->X0.[i,n]=z) N,X0.[i,j])
 X0.[i,j]<-z; randLS i j z i' j' z'

let asNormLS α= 

 let n=Array.init (Array2D.length1 α) (fun n->(α.[n,0]-1,n))|>Map.ofArray
 let g=Array.init (Array2D.length1 α) (fun g->(α.[n.[0],g]-1,g))|>Map.ofArray
 Array2D.init (Array2D.length1 α) (Array2D.length1 α) (fun i j->α.[n.[i],g.[j]])

let randLS α=Seq.unfold(fun g->Some(g,jmLS α g))(Array2D.init α α (fun n g->1+(n+g)%α)) </lang>

The Task

part 1

<lang fsharp> randLS 4 |> Seq.take 10000 |> Seq.map asNormLS |> Seq.countBy id |> Seq.iter(fun n->printf "%A was produced %d times\n\n" (fst n)(snd n)) </lang>

Output:
[[1; 2; 3; 4]
 [2; 3; 4; 1]
 [3; 4; 1; 2]
 [4; 1; 2; 3]] was produced 2920 times

[[1; 2; 3; 4]
 [2; 4; 1; 3]
 [3; 1; 4; 2]
 [4; 3; 2; 1]] was produced 2262 times

[[1; 2; 3; 4]
 [2; 1; 4; 3]
 [3; 4; 2; 1]
 [4; 3; 1; 2]] was produced 2236 times

[[1; 2; 3; 4]
 [2; 1; 4; 3]
 [3; 4; 1; 2]
 [4; 3; 2; 1]] was produced 2582 times
part 2

<lang fsharp> randLS 4 |> Seq.take 10000 |> Seq.map asNormLS |> Seq.countBy id |> Seq.iteri(fun n g->printf "%d(%d) " (n+1) (snd g)); printfn "" </lang>

Output:
1(176) 2(171) 3(174) 4(165) 5(168) 6(182) 7(138) 8(205) 9(165) 10(174) 11(157) 12(187) 13(181) 14(211) 15(184) 16(190) 17(190) 18(192) 19(146) 20(200) 21(162) 22(153) 23(193) 24(156) 25(148) 26(188) 27(186) 28(198) 29(178) 30(217) 31(185) 32(172) 33(223) 34(147) 35(203) 36(167) 37(188) 38(152) 39(165) 40(187) 41(160) 42(199) 43(140) 44(202) 45(186) 46(182) 47(175) 48(161) 49(179) 50(175) 51(201) 52(195) 53(205) 54(183) 55(155) 56(178)
part 3

<lang fsharp> let q=Seq.item 749 (randLS 42) for n in [0..41] do (for g in [0..41] do printf "%3d" q.[n,g]); printfn "" </lang>

Output:
 16  7 41 15 17 40 12  9 10  5 19 29 21 18  8 22  3 36 23 31 11 38 13 30  2 33  6 42 39 14 32 20 28 35 26  1 34 37 27 24  4 25
 38 25 36 32 40 29 35 27  8 26 31 15  9  7 16 11  4  3 12 20 23 33  5 24 41 14 30 34 42 17 39 18 37 22 21 13  1 10  6 19  2 28
  8 34 27 25 21 31  1 23 37 36 26 13 22 24 35 17 10 40 41 30 42  7 15  2 18  3 29 11 32  4 38 39  9  5 16 14 28 12 20 33 19  6
 33 35 13 34 15 24  4 29 41 27  3 17 10 26 39 23 30 32  1 38 16 25 37 14  6 28 19  9 40  5 18  7 42 11 31 20 12 22  2 21  8 36
  2 42 20  1  7 26 11 10 39 41 34 22 40 23 24 29 14 17  5 33 38 30  6 13  3 16 18 19 31 15 28 21 36 37 32 27  8  4 25  9 35 12
 25 33 14 40 28 30 31 24 29  4  8 20 26 38 12 35  2 39 16  6 13 21 18 17  5 41 23  3 36  7 34 22 27  1 10 42 11 19 15 32 37  9
 17 22 35 28 30 18 21  2 15 39  5 40 27 13  1 34 38 37 26 23 41 36  4  3 11  6 20  8  9 10 12 24 31 25  7 29 16 32 42 14 33 19
 14  9 19  7 26 15 10  4 36 25 22 23 39 16  2 40 18  1 38 13 21 37 34 31 35 24 12 27 11  3  5  6 17 20 41 33 32 29  8 30 28 42
  5 27 24 13  2 36 25 30 23  9  6 14 35 15 42 39 16 26 21 34 33 31  3  1 29 12 38 17 37 19 40  4  7  8 22 41 20 28 32 10 18 11
 19 41 28 26  8 10 30 35 18 33 15 27 25 21 29 42 23 12 17  2  5  1 38  6 20  7 34  4 13 36 24 31 14  3 11 32 39 40  9 22 16 37
 41 10  3 19 22  9 27 40  1 29 16 42 33 39 34  7 37 20 11 12  4 18 35  8 28 26 36  5 17 30 25 32  6 15 24 21 13 23 14  2 38 31
 42  3 16 36 33 21 20 14 31 22  9 38 29 19 37 13 28 10 35 18 39 26 25 27  4 30 15 23 41 24 11  1 40  7  5 17  6  2 12  8 34 32
 23 31 34 41 38 33  3 28  4  1 30 25  6  2 20 14 13 24  8 42  7 12 39 32 22 29  5 37 15  9 27 10 35 36 19 40 17 18 16 11 26 21
 37 16 30 11  4 32 42 33 13  6 14  2 15 27 18 31 20 41 39 40  9 24 36  5 10  8  1 26  3 34 22 28 38 19 29 23 21 25 35 12 17  7
  1 19 26 22 16 25 36 39  3 23 41 37 34  6 17 32 40 21 10 27 12  9 31  7 13  4 24 29  8 11  2  5 15 18 35 28 30 20 33 38 42 14
 11 13 23 30 25 41  6 31 14 32 27 36 19 17 10 33 21 15  7  5  8 28 16 35 34 42 40  2 38 39  9 26 20 24 37  4 18  3 22  1 12 29
 24 17 29 38 23 39 32  5 11 15 35 12  8 10 40  1 22 25  2 36 28  4 42 21  9 20  3 31 16 41 13 30 19 34 33 18 27  6  7 37 14 26
 36  4  6 24 12 20  2 34 40 11 32  9 28  8 38 21  5 31 42 17 14 29 19 22 25 15  7 18 30 26  1 13 16 41 23 39 37 33  3 35 10 27
 20 39  2 12 32  7 22  3 17 10 37  6 18 40 27  5 42 35 28  4 24 14 33 29 30 31 26 13 19 23 36 41  1 21  9 11 15  8 34 16 25 38
 35 18 37  6  5 13 29  8 24 19 38 34 12 31 21 10 33  7  3 41 15 42 20 11 27 40 16 14 23  1  4  2 22 32 28  9 25 30 26 39 36 17
 10 32  9 33 39 19 41 38 35 18 28 26 14 30  7  4  1 22 37 21 31 40 27 15 42 34  2 25  5 12 23 36  8  6 17  3 29 24 11 13 20 16
 13 28 39  2 31  8  9 37 21 16 40 19 42 36 41  3 12 14 20 10 17 34  1 33 32 35 25 30 18 38 15 11 24 23  6 26  4  5 29  7 27 22
  7 40 12 39 18  3 16 21 42 17  1 32  5 33 13  6 41  8 29 14 34 35 24 36 38 25 31 28 26 27 20 37 23  2 30 10 22  9 19  4 11 15
  4 21  7 17 35 34 19 25 12 42 11  1 30 28 36 26 32 23 14 29  2 20  8 41 24 27 22 15 10 18 37  9 39 38 13  6  3 16 31 40  5 33
 34 23 42 14 41 27 37  6  9 31  4  5  7  1 25 16 35 30 33 11 19  3 26 12 17 38  8 20 24 13 29 15 32 28 40 22  2 39 18 36 21 10
 30  6 21  9 20 17  5 32 38 13 12 28 16 35 22 36 34 29 40 39 25 15 14 37 33 11  4 41  1  2 19  3 26 27 42  8 10  7 23 31 24 18
  6 38  8 10 42 35 13  1 16 37 21  3 11 34 32 20 29 18 25 22 36  5 30 26 39 23 28 12  2 31  7 19 33 40 14 24  9 41 17 27 15  4
 29 15  1 21 14 11 26 17 30 38 10 33 36 20  4 18 39 16 31  3 35  2 32 28 19 13 42  7 12  8  6 40  5  9 25 37 24 27 41 23 22 34
 21 36 32  8  6 23 15 19  2 14 18  4  3 11  5 28 26 13 34 25 30 17  7 42 16 22 39 40 29 37 33 12 41 10 27 31 35 38 24 20  9  1
 39 20 31 29 19  4 38 16 27 30 24 11  2  3 33 15  8 28 18 37 10 13  9 23 36  1 17 22 25 32 26 35 12 42 34  7 40 14 21  5  6 41
 12 11 17 42  9  2 14  7 22 24 25 31 38 41 15 19 36 33 32 28  1 10 29 40 23 18 37 39  6 21 35 27  3 16  8 30  5 26  4 34 13 20
 18 29 33 16 27 42 40 26  7  8 39 24 41  5 30 38  6  9 13  1 32 22  2 34 12 37 11 10 35 20 14 17 21  4 15 19 23 36 28 25 31  3
 28  2  4 18 11  5 23 20 25 35 42 30 31 14  3  9 24 27 19  7 22  6 12 10  1 32 41 36 21 33 16 34 29 13 39 15 38 17 37 26 40  8
  3 26 11 35 24 37 17 36  6  7 13 41  4 32  9  2 31 34 22 15 29  8 40 18 21  5 27  1 14 16 10 38 25 33 20 12 19 42 39 28 30 23
 31  5 22 27 10  6  8 13 34  2 33  7 32 42 26 12 19  4 15  9 40 16 28 38 37 39 35 24 20 29 17 23 11 14  3 25 41 21 36 18  1 30
 15 24  5 37  3 28  7 22 19 34 20 18 17 12 23  8 25 11 36 16 27 41 10  4 31  2  9 32 33 42 21 14 13 29 38 35 26  1 30  6 39 40
 27 37 25  5 13 16 24 41 28  3  2 10 23  4 14 30 11 38  6 19 26 32 21 20 40  9 33 35 34 22 42  8 18 17 12 36 31 15  1 29  7 39
 26 30 10  3 36 22 33 11  5 20 29 21 13 25 31 37 17  2  9 35 18 27 23 39 14 19 32 16 28  6  8 42  4 12  1 38  7 34 40 15 41 24
 32  8 18 31  1 14 34 12 33 28 17 39 37  9 19 27  7  5 30 24 20 23 11 25 15 36 21  6 22 40 41 16 10 26  4  2 42 35 38  3 29 13
  9 14 40 23 37 38 18 15 20 12 36  8  1 22 28 24 27 42  4 32  6 11 41 19 26 10 13 21  7 25 30 29 34 39  2 16 33 31  5 17  3 35
 22 12 15  4 34  1 39 42 32 40  7 35 20 29 11 25  9  6 24 26 37 19 17 16  8 21 14 38 27 28  3 33 30 31 18  5 36 13 10 41 23  2
 40  1 38 20 29 12 28 18 26 21 23 16 24 37  6 41 15 19 27  8  3 39 22  9  7 17 10 33  4 35 31 25  2 30 36 34 14 11 13 42 32  5

Go

The J & M implementation is based on the C code here which has been heavily optimized following advice and clarification by Nigel Galloway (see Talk page) on the requirements of this task.

Part 4 is taking about 6.5 seconds on my Celeron @1.6 GHz but will be much faster on a more modern machine. Being able to compute random, uniformly distributed, Latin squares of order 256 reasonably quickly is interesting from a secure communications or cryptographic standpoint as the symbols of such a square can represent the 256 characters of the various extended ASCII encodings. <lang go>package main

import ( 

   "fmt"
   "math/rand"
   "time"

)

type ( 

   vector []int
   matrix []vector
   cube   []matrix

)

func toReduced(m matrix) matrix { 

   n := len(m)
   r := make(matrix, n)
   for i := 0; i < n; i++ {
       r[i] = make(vector, n)
       copy(r[i], m[i])
   }
   for j := 0; j < n-1; j++ {
       if r[0][j] != j {
           for k := j + 1; k < n; k++ {
               if r[0][k] == j {
                   for i := 0; i < n; i++ {
                       r[i][j], r[i][k] = r[i][k], r[i][j]
                   }
                   break
               }
           }
       }
   }
   for i := 1; i < n-1; i++ {
       if r[i][0] != i {
           for k := i + 1; k < n; k++ {
               if r[k][0] == i {
                   for j := 0; j < n; j++ {
                       r[i][j], r[k][j] = r[k][j], r[i][j]
                   }
                   break
               }
           }
       }
   }
   return r

}

// 'm' is assumed to be 0 based func printMatrix(m matrix) { 

   n := len(m)
   for i := 0; i < n; i++ {
       for j := 0; j < n; j++ {
           fmt.Printf("%2d ", m[i][j]+1) // back to 1 based
       }
       fmt.Println()
   }
   fmt.Println()

}

// converts 4 x 4 matrix to 'flat' array func asArray16(m matrix) [16]int { 

   var a [16]int
   k := 0
   for i := 0; i < 4; i++ {
       for j := 0; j < 4; j++ {
           a[k] = m[i][j]
           k++
       }
   }
   return a

}

// converts 5 x 5 matrix to 'flat' array func asArray25(m matrix) [25]int { 

   var a [25]int
   k := 0
   for i := 0; i < 5; i++ {
       for j := 0; j < 5; j++ {
           a[k] = m[i][j]
           k++
       }
   }
   return a

}

// 'a' is assumed to be 0 based func printArray16(a [16]int) { 

   for i := 0; i < 4; i++ {
       for j := 0; j < 4; j++ {
           k := i*4 + j
           fmt.Printf("%2d ", a[k]+1) // back to 1 based
       }
       fmt.Println()
   }
   fmt.Println()

}

func shuffleCube(c cube) { 

   n := len(c[0])
   proper := true
   var rx, ry, rz int
   for {
       rx = rand.Intn(n)
       ry = rand.Intn(n)
       rz = rand.Intn(n)
       if c[rx][ry][rz] == 0 {
           break
       }
   }
   for {
       var ox, oy, oz int
       for ; ox < n; ox++ {
           if c[ox][ry][rz] == 1 {
               break
           }
       }
       if !proper && rand.Intn(2) == 0 {
           for ox++; ox < n; ox++ {
               if c[ox][ry][rz] == 1 {
                   break
               }
           }
       }

       for ; oy < n; oy++ {
           if c[rx][oy][rz] == 1 {
               break
           }
       }
       if !proper && rand.Intn(2) == 0 {
           for oy++; oy < n; oy++ {
               if c[rx][oy][rz] == 1 {
                   break
               }
           }
       }

       for ; oz < n; oz++ {
           if c[rx][ry][oz] == 1 {
               break
           }
       }
       if !proper && rand.Intn(2) == 0 {
           for oz++; oz < n; oz++ {
               if c[rx][ry][oz] == 1 {
                   break
               }
           }
       }

       c[rx][ry][rz]++
       c[rx][oy][oz]++
       c[ox][ry][oz]++
       c[ox][oy][rz]++

       c[rx][ry][oz]--
       c[rx][oy][rz]--
       c[ox][ry][rz]--
       c[ox][oy][oz]--

       if c[ox][oy][oz] < 0 {
           rx, ry, rz = ox, oy, oz
           proper = false
       } else {
           proper = true
           break
       }
   }

}

func toMatrix(c cube) matrix { 

   n := len(c[0])
   m := make(matrix, n)
   for i := 0; i < n; i++ {
       m[i] = make(vector, n)
   }
   for i := 0; i < n; i++ {
       for j := 0; j < n; j++ {
           for k := 0; k < n; k++ {
               if c[i][j][k] != 0 {
                   m[i][j] = k
                   break
               }
           }
       }
   }
   return m

}

// 'from' matrix is assumed to be 1 based func makeCube(from matrix, n int) cube { 

   c := make(cube, n)
   for i := 0; i < n; i++ {
       c[i] = make(matrix, n)
       for j := 0; j < n; j++ {
           c[i][j] = make(vector, n)
           var k int
           if from == nil {
               k = (i + j) % n
           } else {
               k = from[i][j] - 1
           }
           c[i][j][k] = 1
       }
   }
   return c

}

func main() { 

   rand.Seed(time.Now().UnixNano())

   // part 1
   fmt.Println("PART 1: 10,000 latin Squares of order 4 in reduced form:\n")
   from := matrix{{1, 2, 3, 4}, {2, 1, 4, 3}, {3, 4, 1, 2}, {4, 3, 2, 1}}
   freqs4 := make(map[[16]int]int, 10000)
   c := makeCube(from, 4)
   for i := 1; i <= 10000; i++ {
       shuffleCube(c)
       m := toMatrix(c)
       rm := toReduced(m)
       a16 := asArray16(rm)
       freqs4[a16]++
   }
   for a, freq := range freqs4 {
       printArray16(a)
       fmt.Printf("Occurs %d times\n\n", freq)
   }

   // part 2
   fmt.Println("\nPART 2: 10,000 latin squares of order 5 in reduced form:")
   from = matrix{{1, 2, 3, 4, 5}, {2, 3, 4, 5, 1}, {3, 4, 5, 1, 2},
       {4, 5, 1, 2, 3}, {5, 1, 2, 3, 4}}
   freqs5 := make(map[[25]int]int, 10000)
   c = makeCube(from, 5)
   for i := 1; i <= 10000; i++ {
       shuffleCube(c)
       m := toMatrix(c)
       rm := toReduced(m)
       a25 := asArray25(rm)
       freqs5[a25]++
   }
   count := 0
   for _, freq := range freqs5 {
       count++
       if count > 1 {
           fmt.Print(", ")
       }
       if (count-1)%8 == 0 {
           fmt.Println()
       }
       fmt.Printf("%2d(%3d)", count, freq)
   }
   fmt.Println("\n")

   // part 3
   fmt.Println("\nPART 3: 750 latin squares of order 42, showing the last one:\n")
   var m42 matrix
   c = makeCube(nil, 42)
   for i := 1; i <= 750; i++ {
       shuffleCube(c)
       if i == 750 {
           m42 = toMatrix(c)
       }
   }
   printMatrix(m42)

   // part 4
   fmt.Println("\nPART 4: 1000 latin squares of order 256:\n")
   start := time.Now()
   c = makeCube(nil, 256)
   for i := 1; i <= 1000; i++ {
       shuffleCube(c)
   }
   elapsed := time.Since(start)
   fmt.Printf("Generated in %s\n", elapsed)

}</lang>

Output:

Sample run: 

PART 1: 10,000 latin Squares of order 4 in reduced form:

 1  2  3  4 
 2  1  4  3 
 3  4  2  1 
 4  3  1  2 

Occurs 2550 times

 1  2  3  4 
 2  4  1  3 
 3  1  4  2 
 4  3  2  1 

Occurs 2430 times

 1  2  3  4 
 2  1  4  3 
 3  4  1  2 
 4  3  2  1 

Occurs 2494 times

 1  2  3  4 
 2  3  4  1 
 3  4  1  2 
 4  1  2  3 

Occurs 2526 times


PART 2: 10,000 latin squares of order 5 in reduced form:

 1(165),  2(173),  3(167),  4(204),  5(173),  6(165),  7(215),  8(218), 
 9(168), 10(157), 11(205), 12(152), 13(187), 14(173), 15(215), 16(185), 
17(179), 18(176), 19(179), 20(160), 21(150), 22(166), 23(191), 24(181), 
25(179), 26(192), 27(187), 28(186), 29(176), 30(196), 31(141), 32(187), 
33(165), 34(189), 35(147), 36(175), 37(172), 38(162), 39(180), 40(172), 
41(189), 42(159), 43(197), 44(158), 45(178), 46(179), 47(193), 48(175), 
49(207), 50(174), 51(181), 52(179), 53(193), 54(171), 55(153), 56(204)


PART 3: 750 latin squares of order 42, showing the last one:

29  2 17 41 34 30  8 33 39  7 20 27 12  6 31 14 40 35 25  9 10 32 19 16 24 42  3 26  5 23  1 28  4 13 38 18 21 37 22 15 36 11 
17 15 11 31  9 38 26 10  1 28 37  8 34 41 21 22 12  5 35 36 13 20 29 42 18  3 19 24 39 32 27 23 16 25 33  4 40  6  2 30  7 14 
36 42 35 39 15 34 37 18 32 25 22 31  4 17  3 19 13 11  8 23 12 24 28 27 16  1  6  9 29 40  7  5  2 14 30 26 41 10 21 33 38 20 
21 13 16 42  3 32  2 26 27 17 15 11 25 37 29  6 19 10 12  7 31 18 36  9 39 41 30 40 35 33 22  1 28 38 24  8 34 23  4 20 14  5 
22 39 13  7 38  9 34 41 37 36 35  6 21 26 17 16  4 30 40 20  8 15 25 19 32  2 11 28 23 24 31 10 42  3 27 12 33 14  1 29  5 18 
33 36 34  3 13  4  7 14  2 29  6 12 31 23 26 17  8 20 32 21 19 41 37  5 38 30 25 11 24 35 42 27 18 16 39 15 10 22 28  1  9 40 
14 31  7 22 39 23 32 34 16 33 24  4 40 42 12 25 35 26 18 28 11  3 15 21 20  9 13 19  1 10  2 41 29  6 17 30  5 38 37  8 27 36 
 9  3  6 30 19 39 14 16  4 15 29 28 23 24 32 10 18 41 37 38 40 34  8 25  2 22 31  5 17 26 36 33 13 21 12 35  7 20 11 27 42  1 
 2 18 28  5  6  7 40 35  3 20  8 34 42 39 37 33 26 23 22 13 14  4 12 15 17 25 36 31 16 29 38 19 32 41  1 27 24 11 30  9 10 21 
27 34 19 15 33 22  5 36  9 30 14  1 24  8 38 42 41 39  7 40  4 37 11 23 29 26 18 12  3 21 35 16 20 10 31 25 17 28  6 32  2 13 
41 16  1 35 22 13 20 29  6 38  5 24 19 10 25 27 17 18 11 32  9  7  2 36  4 34 40 21 33 12  8 30 15 42 37 23 14 26  3 39 31 28 
 7  1 15 16 27 31 18 24 20  8 36 38 10 34  9  4 42 29  2  3 26 39  5 22 41 21 37 30 14 11 33 35 25 23 40 28 13 19 17  6 32 12 
 1 10 20 32 23  5 30 12  8  9 21 36 15 14 18 37 33 31 26 39 41 16  6 24 22 35 29 42 27 28  3 38 11  2  7 34  4 40 19 17 13 25 
 6 32 42 11 20 40 27 25 41 22 17 16 26 29 15  7 23 36 39 34 28 13 18  3 10 37  8 14  2 31  4 24  5 19  9 21 38  1 33 12 30 35 
35 40 30 19 21 12 17  4 22 27  3 20 11  9  8 23 24 42 14 10 39 28 26 29 33 13 41 16 34 25 32 37  7 18  5  6 15  2 36 38  1 31 
15 26 40  1 28 20  9 21  7  5 13 18 30 22 10  8  3 25  6  2 17 36 38 31 14 19 35 23 12 27 11 39 24  4 41 32 29 34 42 16 37 33 
 3  6 26 12 32  1 13  8 42 37 25  7  9 16 35  5 29 21 24 27 34 17 14  2 15 11 28 33 20 38 18 22 39 40 23 10 31 30 41 36 19  4 
31 38 36 21 16 26 28 30 15  3 32 41 18  1  6 29  9 17  5 35  7 40 27 37 13 20 23 22 11 19 12 42 34  8 10 14 25 39 24  4 33  2 
40  4 22 38 35 11 21 17 31  1 28 19 37  2 42 24 14 12 13 30 33 25 34 32 27 36 39  3  9 15 10 18  8  5  6 41 26 16 29  7 20 23 
 5 17 39  4 26 14 31 37 35 11 38  3  1 30 19 36 20 33 15 16 21 29  9  6 25 27  2 13 41 34 24 12 10 32 22  7 28 18 40 42 23  8 
 8 29 24 26 31 21 39 23 11 14 19 10 20 15  7 35 32 38  1 12 25 22 16  4  6 40 42 41 18 30 28  2 17 36  3 13 37 33 27  5 34  9 
11 25 14 17 18 24 19 32 33 31  7 26  2 21 20 30 15 27 23 41 29 35 39 28 34 12 10  4  8 42  5 13 37  9 16 40  1 36 38  3  6 22 
26 21 18 25 29 15  1 13 19  2 34 23 38 27 41  3 10 22 17  4 16 11 42 12  8  6  5 35 30 39 37 14  9 24 36 33 20  7 31 28 40 32 
25 27 12 33 17 35 24  9 28 10 42 21  8 13  2 15 34 16  3 18  5 31 41  7 23  4  1  6 22 14 19 36 40 37 26 38 30 32 20 11 39 29 
23 19 25  9 30 37 38 40 14 41 31 17  7  4 16 11  1  6 33  5 24  2  3  8 21 29 34 32 28 22 15 20 12 35 18 36 39 27 10 13 26 42 
34  9 10 13  2  6 22 31 26 40  1 14 41  3 11 12 37 32 27 29 35 19 30 33 28 38 21 25  7  5 16  8 36 15 20 42 23 17 39 18  4 24 
20 11 37 28 41  8 10 15 36 12 26 33 39 32 13  1 25  9 42 19  3  6 24 14  5 23  7 27 38  2 30  4 22 34 35 31 18 29 16 40 21 17 
28 30 21 23 24 29  3  1 10  6 33  2 27 40 14 34 31 15 19 37 18  9  4 13 35  8 12 20 36 16 17 32 41  7 25 39 42  5 26 22 11 38 
32 12  8 40 11 16 23 28 18 42 41 30  3 38 33  2 22 19  4 25 37  1 31 20 36  5  9  7 13 17 14  6 27 39 34 24 35 21 15 26 29 10 
18 37 41 10 36 28 11 42 13 34  2 35  5  7 22 40 39  3 30  1 38 27 20 17 19 33 26 15 25  6 21 29 23 31  4  9 32  8 12 14 24 16 
39 24 29 37 25 19 33 27 17 16 10 40 36 12 30 41 11  4 34 15  2  5 32  1 31 14 38 18 42  3  9  7  6 20 21 22  8 13 23 35 28 26 
19 14  5  8 40  3 29  6 21 26 23 15 16 33 28 31 38 13  9 17 27 12 10 11  7 24 20  1  4 41 39 25 30 22 32  2 36 42 35 34 18 37 
37  7 32 34  8 36 41  2 12 24 16 39 33 31  4 13  6 28 38 22 20 42 40 18  9 10 14 29 26  1 23 15 21 27 19 17 11  3  5 25 35 30 
 4 41 27  2 42 17 15 38 30 35 12 25 13 28 39 20  5  1 16 33 36 23  7 40 37 32 24 10 31  8  6 21 14 26 29 11  3  9 18 19 22 34 
38 35 23 36  4 10 12 11  5 21 27 32 17 25 24 18 28 40 20  6 42 14 22 30 26 39 33  8 37  7 13 34  1 29 15 19  2 41  9 31 16  3 
30 33 31 24 12 41 36 19 23 32  4 37 29 11 34 39 16 14 21 42  6 26  1 38  3 17 22  2 40 18 20  9 35 28 13  5 27 15 25 10  8  7 
42 28  3 14  1 25 16 22 34 23 39  9 35  5 40 26 36  7 10 31 32 21 13 41 30 18  4 38  6 37 29 17 33 12 11 20 19 24  8  2 15 27 
16  5 38  6 10 27  4  3 40 18 11 13 22 35  1 21  2 34 36  8 23 30 17 39 42  7 15 37 32 20 26 31 19 33 28 29  9 25 14 24 12 41 
24 23 33 18 14  2 25 39 29 19  9  5 28 20 27 38  7  8 31 11 15 10 35 34 12 16 32 17 21 36 40  3 26 30 42  1 22  4 13 37 41  6 
12 20  2 29  5 33 42  7 24  4 18 22 14 19 36  9 27 37 28 26 30 38 23 10 11 31 17 34 15 13 41 40  3  1  8 16  6 35 32 21 25 39 
13  8  9 27 37 42  6 20 25 39 40 29 32 18  5 28 30 24 41 14 22 33 21 35  1 15 16 36 10  4 34 26 38 11  2  3 12 31  7 23 17 19 
10 22  4 20  7 18 35  5 38 13 30 42  6 36 23 32 21  2 29 24  1  8 33 26 40 28 27 39 19  9 25 11 31 17 14 37 16 12 34 41  3 15 


PART 4: 1000 latin squares of order 256:

Generated in 6.581088256s
Retrieved from "https://rosettacode.org/wiki/Latin_Squares_in_reduced_form/Randomizing_using_Jacobson_and_Matthews’_Technique?oldid=292155"