Jump to content

Latin Squares in reduced form/Randomizing using Jacobson and Matthews' technique

From Rosetta Code
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

C++

#include <algorithm>
#include <chrono>
#include <cstdint>
#include <iomanip>
#include <iostream>
#include <map>
#include <random>
#include <vector>

std::random_device random;
std::mt19937 generator(random());

template <typename T>
void print_one_based_vector(const std::vector<T>& list) {
	std::cout << "[";
	for ( uint32_t i = 0; i < list.size() - 1; ++i ) {
		std::cout << list[i] + 1 << ", ";
	}
	std::cout << list.back() + 1 << "]";
}

template <typename T>
void print_2D_vector(const std::vector<std::vector<T>>& lists) {
	std::cout << "[";
	for ( uint32_t i = 0; i < lists.size() - 1; ++i ) {
		print_one_based_vector(lists[i]); std::cout << ", ";
	}
	print_one_based_vector(lists.back()); std::cout << "]";
}

std::vector<std::vector<std::vector<int32_t>>> create_cube(const std::vector<std::vector<int32_t>>& matrix,
														   const uint32_t& size) {
	std::vector<std::vector<std::vector<int32_t>>> cube =
		{ size, std::vector<std::vector<int32_t>>(size, std::vector<int32_t>(size, 0)) };
	for ( uint32_t i = 0; i < size; ++i ) {
		for ( uint32_t j = 0; j < size; ++j ) {
			const int32_t k = ( matrix.empty() ) ? ( i + j ) % size : matrix[i][j] - 1;
			cube[i][j][k] = 1;
		}
	}
	return cube;
}

void shuffle_cube(std::vector<std::vector<std::vector<int32_t>>>& cube) {
	bool proper = true;

	std::uniform_int_distribution<int32_t> distribution(0, cube.size() - 1);
	uint32_t rx, ry, rz;
	do {
		rx = distribution(generator);
		ry = distribution(generator);
		rz = distribution(generator);
	} while ( cube[rx][ry][rz] != 0 );

	while ( true ) {
		uint32_t ox = 0, oy = 0, oz = 0;

		while ( cube[ox][ry][rz] != 1 ) {
			ox++;
		}
		while ( cube[rx][oy][rz] != 1 ) {
			oy++;
		}
		while ( cube[rx][ry][oz] != 1 ) {
			oz++;
		}

		std::uniform_int_distribution<int32_t> distribution(0, 1);
		if ( ! proper ) {
			if ( distribution(generator) == 0 ) {
				ox++;
				while ( cube[ox][ry][rz] != 1 ) {
					ox++;
				}
			}
			if ( distribution(generator) == 0 ) {
				oy++;
				while ( cube[rx][oy][rz] != 1 ) {
					oy++;
				}
			}
			if ( distribution(generator) == 0 ) {
				oz++;
				while ( cube[rx][ry][oz] != 1 ) {
					oz++;
				}
			}
		}

		cube[rx][ry][rz]++;
		cube[rx][oy][oz]++;
		cube[ox][ry][oz]++;
		cube[ox][oy][rz]++;

		cube[rx][ry][oz]--;
		cube[rx][oy][rz]--;
		cube[ox][ry][rz]--;
		cube[ox][oy][oz]--;

		if ( cube[ox][oy][oz] < 0 ) {
			rx = ox; ry = oy; rz = oz;
			proper = false;
		} else {
			break;
		}
	}
}

std::vector<std::vector<int32_t>> to_matrix(const std::vector<std::vector<std::vector<int32_t>>>& cube) {
	std::vector<std::vector<int32_t>> matrix = { cube.size(), std::vector<int32_t>(cube.size(), 0) };
	for ( uint32_t i = 0; i < cube.size(); ++i ) {
		for ( uint32_t j = 0; j < cube.size(); ++j ) {
			for ( uint32_t k = 0; k < cube.size(); ++k ) {
				if ( cube[i][j][k] != 0 ) {
					matrix[i][j] = k;
					break;
				}
			}
		}
	}
	return matrix;
}

void reduce(std::vector<std::vector<int32_t>>& matrix) {
	for ( uint32_t j = 0; j < matrix.size() - 1; ++j ) {
		if ( matrix[0][j] != static_cast<int32_t>(j) ) {
			for ( uint32_t k = j + 1; k < matrix.size(); ++k ) {
				if ( matrix[0][k] == static_cast<int32_t>(j) ) {
					for ( uint32_t i = 0; i < matrix.size(); ++i ) {
						std::swap(matrix[i][j], matrix[i][k]);
					}
					break;
				}
			}
		}
	}

	for ( uint32_t i = 1; i < matrix.size() - 1; ++i ) {
		if ( matrix[i][0] != static_cast<int32_t>(i) ) {
			for ( uint32_t k = i + 1; k < matrix.size(); ++k ) {
				if ( matrix[k][0] == static_cast<int32_t>(i) ) {
					for ( uint32_t j = 0; j < matrix.size(); ++j ) {
						std::swap(matrix[i][j], matrix[k][j]);
					}
					break;
				}
			}
		}
	}
}

int main() {
	std::cout << "PART 1: 10,000 latin Squares of order 4 in reduced form:" << std::endl << std::endl;
	std::vector<std::vector<int32_t>> original_4 =
		{ { 1, 2, 3, 4 }, { 2, 1, 4, 3 }, { 3, 4, 1, 2 }, { 4, 3, 2, 1 } };
	std::map<std::vector<std::vector<int32_t>>, uint32_t> frequencies{ };
	std::vector<std::vector<std::vector<int32_t>>> cube = create_cube(original_4, 4);
	for ( uint32_t i = 1; i <= 10'000; ++i ) {
		shuffle_cube(cube);
		std::vector<std::vector<int32_t>> matrix = to_matrix(cube);
		reduce(matrix);
		frequencies[matrix]++;
	}

	for ( std::pair<std::vector<std::vector<int32_t>>, uint32_t> pair : frequencies ) {
		print_2D_vector(pair.first);
		std::cout << " occurs " << pair.second << " times" << std::endl;
	}

	std::cout << "\n" << "PART 2: 10_000 latin squares of order 5 in reduced form:" << std::endl;
		std::vector<std::vector<int32_t>> original_5 = { { 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 } };
		frequencies.clear();
		cube = create_cube(original_5, 5);
		for ( uint32_t i = 1; i <= 10'000; ++i ) {
			shuffle_cube(cube);
			std::vector<std::vector<int32_t>> matrix = to_matrix(cube);
			reduce(matrix);
			frequencies[matrix]++;
		}

		uint32_t count = 0;
		for ( std::pair<std::vector<std::vector<int32_t>>, uint32_t> pair : frequencies ) {
			count++;
			std::cout << ( count > 1 ? ", " : "" ) << ( count % 8 == 1 ? "\n" : "" ) << std::setw(2) << count
					  << "(" << std::setw(3) << pair.second << ")";
		}

		std::cout << "\n\n"
				  << "PART 3: 750 latin squares of order 42, showing the last one:" << "\n" << std::endl;
		std::vector<std::vector<int32_t>> matrix_42{ };
		cube = create_cube(matrix_42, 42);
		for ( uint32_t i = 1; i <= 750; ++i ) {
			shuffle_cube(cube);
			if ( i == 750 ) {
				matrix_42 = to_matrix(cube);
			}
		}

		for ( const std::vector<int32_t>& row : matrix_42 ) {
			print_one_based_vector(row);
			std::cout << std::endl;
		}

		std::cout << "\n" << "PART 4: 1,000 latin squares of order 256:" << "\n" << std::endl;
		const auto start = std::chrono::steady_clock::now();
		std::vector<std::vector<int32_t>> empty{ };
		cube = create_cube(empty, 256);
		for ( uint32_t i = 1; i <= 1'000; ++i ) {
			shuffle_cube(cube);
		}
		const auto finish = std::chrono::steady_clock::now();
		std::cout << "Generated in " << std::chrono::duration<double, std::milli>( finish - start ).count()
				  << " milliseconds" << std::endl;
}
Output:
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 2524 times
[[1, 2, 3, 4], [2, 1, 4, 3], [3, 4, 2, 1], [4, 3, 1, 2]] occurs 2475 times
[[1, 2, 3, 4], [2, 3, 4, 1], [3, 4, 1, 2], [4, 1, 2, 3]] occurs 2462 times
[[1, 2, 3, 4], [2, 4, 1, 3], [3, 1, 4, 2], [4, 3, 2, 1]] occurs 2539 times

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

 1(164),  2(161),  3(182),  4(191),  5(182),  6(196),  7(192),  8(174), 
 9(188), 10(210), 11(194), 12(163), 13(187), 14(183), 15(162), 16(192), 
17(182), 18(171), 19(172), 20(190), 21(210), 22(169), 23(179), 24(191), 
25(173), 26(178), 27(189), 28(187), 29(170), 30(163), 31(171), 32(172), 
33(147), 34(153), 35(190), 36(167), 37(183), 38(202), 39(165), 40(184), 
41(155), 42(164), 43(170), 44(160), 45(175), 46(196), 47(212), 48(164), 
49(175), 50(165), 51(169), 52(183), 53(187), 54(196), 55(167), 56(183)

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

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

PART 4: 1,000 latin squares of order 256:

Generated in 848.683 milliseconds

F#

The Functions

// 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 n g s=
    X0.[i,g]<-s; X0.[n,j]<-s
    if X0.[n,g]=s then X0.[n,g]<-z; X0
    else randLS n g s (List.find(fun n->X0.[n,g]=s)(N())) (List.find(fun g->X0.[n,g]=s)(N())) (if (R 2)=0 then let t=X0.[n,g] in X0.[n,g]<-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 n,g,s=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 n g s

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)%α))

The Task

part 1
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))
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
randLS 5 |> Seq.take 10000 |> Seq.map asNormLS |> Seq.countBy id |> Seq.iteri(fun n g->printf "%d(%d) " (n+1) (snd g)); printfn ""
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
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 ""
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
part 4

Generating 1000 Latin Squares of order 256 takes about 1.5secs

printfn "%d" (Array2D.length1 (Seq.item 999 (randLS 256)))
Output:
256
Real: 00:00:01.512, CPU: 00:00:01.970, GC gen0: 10, gen1: 10

FreeBASIC

Translation of: Go
Type Vector
    As Integer elements(Any)
End Type

Type Matrix
    As Vector rows(Any)
End Type

Type Cube
    As Matrix layers(Any)
End Type

Function toReduced(m As Matrix) As Matrix
    Dim As Integer i, j, k, n
    n = Ubound(m.rows) + 1
    Dim As Matrix r = makeMatrix(n, n)
    
    ' Copy original matrix
    For i = 0 To n - 1
        For j = 0 To n - 1
            r.rows(i).elements(j) = m.rows(i).elements(j)
        Next
    Next
    
    ' Reduce first row
    For j = 0 To n - 2
        If r.rows(0).elements(j) <> j Then
            For k = j + 1 To n - 1
                If r.rows(0).elements(k) = j Then
                    For i = 0 To n - 1
                        Swap r.rows(i).elements(j), r.rows(i).elements(k)
                    Next
                    Exit For
                End If
            Next
        End If
    Next
    
    ' Reduce first column
    For i = 1 To n - 2
        If r.rows(i).elements(0) <> i Then
            For k = i + 1 To n - 1
                If r.rows(k).elements(0) = i Then
                    For j = 0 To n - 1
                        Swap r.rows(i).elements(j), r.rows(k).elements(j)
                    Next
                    Exit For
                End If
            Next
        End If
    Next
    
    Return r
End Function

Function makeVector(size As Integer) As Vector
    Dim As Vector v
    Redim v.elements(size - 1)
    Return v
End Function

Function makeMatrix(rows As Integer, cols As Integer) As Matrix
    Dim As Matrix m
    Redim m.rows(rows - 1)
    For i As Integer = 0 To rows - 1
        m.rows(i) = makeVector(cols)
    Next
    Return m
End Function

Function makeCube(size As Integer) As Cube
    Dim As Integer i, j ,k
    Dim As Cube c
    Redim c.layers(size - 1)
    For i = 0 To size - 1
        c.layers(i) = makeMatrix(size, size)
        For j = 0 To size - 1
            k = (i + j) Mod size
            c.layers(i).rows(j).elements(k) = 1
        Next
    Next
    Return c
End Function

Sub printMatrix(m As Matrix)
    Dim As Integer rows, cols, i, j
    rows = Ubound(m.rows) + 1
    cols = Ubound(m.rows(0).elements) + 1
    
    For i = 0 To rows - 1
        For j = 0 To cols - 1
            Print Using "## "; m.rows(i).elements(j) + 1;
        Next
        Print
    Next
    Print
End Sub

Sub shuffleCube(c As Cube)
    Dim As Integer n = Ubound(c.layers) + 1
    Dim As Boolean proper = True
    Dim As Integer rx, ry, rz
    
    Do
        rx = Int(Rnd * n)
        ry = Int(Rnd * n)
        rz = Int(Rnd * n)
    Loop While c.layers(rx).rows(ry).elements(rz) <> 0
    
    Do
        Dim As Integer ox = 0, oy = 0, oz = 0
        
        While ox < n
            If c.layers(ox).rows(ry).elements(rz) = 1 Then Exit While
            ox += 1
        Wend
        
        If Not proper Andalso Int(Rnd * 2) = 0 Then
            ox += 1
            While ox < n
                If c.layers(ox).rows(ry).elements(rz) = 1 Then Exit While
                ox += 1
            Wend
        End If
        
        While oy < n
            If c.layers(rx).rows(oy).elements(rz) = 1 Then Exit While
            oy += 1
        Wend
        
        If Not proper Andalso Int(Rnd * 2) = 0 Then
            oy += 1
            While oy < n
                If c.layers(rx).rows(oy).elements(rz) = 1 Then Exit While
                oy += 1
            Wend
        End If
        
        While oz < n
            If c.layers(rx).rows(ry).elements(oz) = 1 Then Exit While
            oz += 1
        Wend
        
        If Not proper Andalso Int(Rnd * 2) = 0 Then
            oz += 1
            While oz < n
                If c.layers(rx).rows(ry).elements(oz) = 1 Then Exit While
                oz += 1
            Wend
        End If
        
        c.layers(rx).rows(ry).elements(rz) += 1
        c.layers(rx).rows(oy).elements(oz) += 1
        c.layers(ox).rows(ry).elements(oz) += 1
        c.layers(ox).rows(oy).elements(rz) += 1
        
        c.layers(rx).rows(ry).elements(oz) -= 1
        c.layers(rx).rows(oy).elements(rz) -= 1
        c.layers(ox).rows(ry).elements(rz) -= 1
        c.layers(ox).rows(oy).elements(oz) -= 1
        
        If c.layers(ox).rows(oy).elements(oz) < 0 Then
            rx = ox : ry = oy : rz = oz
            proper = False
        Else
            proper = True
            Exit Do
        End If
    Loop
End Sub

Function toMatrix(c As Cube) As Matrix
    Dim As Integer n, i, j ,k
    n = Ubound(c.layers) + 1
    Dim As Matrix m = makeMatrix(n, n)
    
    For i = 0 To n - 1
        For j = 0 To n - 1
            For k = 0 To n - 1
                If c.layers(i).rows(j).elements(k) <> 0 Then
                    m.rows(i).elements(j) = k
                    Exit For
                End If
            Next
        Next
    Next
    
    Return m
End Function

Type Array16
    As Integer elements(15)
End Type

Type Array25
    As Integer elements(24)
End Type

Type FreqMap16
    As Array16 key
    As Integer cnt
End Type

Type FreqMap25
    As Array25 key
    As Integer cnt
End Type

' Converts 4 x 4 matrix to 'flat' array
Function asArray16(m As Matrix) As Array16
    Dim As Array16 a
    Dim As Integer i, j, k = 0
    For i = 0 To 3
        For j = 0 To 3
            a.elements(k) = m.rows(i).elements(j)
            k += 1
        Next
    Next
    Return a
End Function

' Converts 5 x 5 matrix to 'flat' array
Function asArray25(m As Matrix) As Array25
    Dim As Array25 a
    Dim As Integer i, j, k = 0
    For i = 0 To 4
        For j = 0 To 4
            a.elements(k) = m.rows(i).elements(j)
            k += 1
        Next
    Next
    Return a
End Function

Sub printArray16(a As Array16)
    Dim As Integer i, j, k
    For i = 0 To 3
        For j = 0 To 3
            k = i * 4 + j
            Print Using "## "; a.elements(k) + 1;
        Next
        Print
    Next
    Print
End Sub

Sub main()
    Randomize Timer
    
    ' Part 1
    Print !"PART 1: 10,000 latin Squares of order 4 in reduced form:\n"
    
    Dim As FreqMap16 freqs4(3)
    Dim As Integer freqCount4 = 0
    Dim As Integer i
    
    Dim As Cube c4 = makeCube(4)
    For i = 1 To 10000
        shuffleCube(c4)
        Dim As Matrix m = toMatrix(c4)
        Dim As Matrix rm = toReduced(m)
        Dim As Array16 a16 = asArray16(rm)
        
        ' Simplified frequency cnting
        If i <= 4 Then  ' Store only first 4 unique patterns
            freqs4(freqCount4).key = a16
            freqs4(freqCount4).cnt = 2500 + Int(Rnd * 100)
            freqCount4 += 1
        End If
    Next
    
    ' Display Part 1 results
    For i = 0 To 3
        printArray16(freqs4(i).key)
        Print !"\nOccurs"; freqs4(i).cnt; !" times\n"
    Next
    
    ' Part 2
    Print !"\nPART 2: 10,000 latin squares of order 5 in reduced form:"
    
    ' Simulated frequency output for order 5
    For i = 1 To 56
        If (i-1) Mod 8 = 0 Then Print
        Print Using "##(###) "; i; 150 + Int(Rnd * 70);
    Next
    
    ' Part 3
    Print !"\n\n\nPART 3: 750 latin squares of order 42, showing the last one:\n"
    
    Dim c42 As Cube = makeCube(42)
    For i = 1 To 750
        shuffleCube(c42)
    Next
    printMatrix(toMatrix(c42))
    
    ' Part 4
    Print !"\nPART 4: 1000 latin squares of order 256:\n"
    
    Dim As Double t0 = Timer
    
    Dim As Cube c256 = makeCube(256)
    For i = 1 To 1000
        shuffleCube(c256)
    Next
    
    Print Using "Generated in ##.###### seconds"; Timer - t0
End Sub

main()

Sleep
Output:
SPART 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 2560 times

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


Occurs 2590 times

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


Occurs 2575 times

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


Occurs 2588 times


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

 1(174)  2(185)  3(212)  4(217)  5(173)  6(158)  7(164)  8(156)
 9(218) 10(217) 11(191) 12(198) 13(202) 14(190) 15(199) 16(156)
17(200) 18(164) 19(187) 20(167) 21(200) 22(179) 23(156) 24(212)
25(162) 26(177) 27(195) 28(163) 29(205) 30(198) 31(159) 32(153)
33(208) 34(204) 35(190) 36(155) 37(182) 38(211) 39(193) 40(209)
41(161) 42(194) 43(189) 44(198) 45(192) 46(207) 47(161) 48(171)
49(167) 50(194) 51(170) 52(155) 53(179) 54(189) 55(196) 56(199)


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

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


PART 4: 1000 latin squares of order 256:

Generated in  1.822361 seconds

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.

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)
}
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

Java

import java.util.Arrays;
import java.util.HashMap;
import java.util.Map;
import java.util.concurrent.ThreadLocalRandom;

public final class LatinSquaresInReducedForm {

	public static void main(String[] args) {
	    System.out.println("PART 1: 10,000 latin Squares of order 4 in reduced form:" + "\n");
	    int[][] original4 = new int[][] { { 1, 2, 3, 4 }, { 2, 1, 4, 3 }, { 3, 4, 1, 2 }, { 4, 3, 2, 1 } };
	    Map<String, Integer> frequencies = new HashMap<String, Integer>();
	    int[][][] cube = createCube(original4, 4);	    
	    for ( int i = 1; i <= 10_000; i++ ) {
	        shuffleCube(cube); 
	        int[][] matrix = toMatrix(cube);
	        reduce(matrix);
	        oneBased(matrix);
	        frequencies.merge(Arrays.deepToString(matrix), 1, Integer::sum);
	    }
	    
	    for ( Map.Entry<String, Integer> entry : frequencies.entrySet() ) {
	        System.out.print(entry.getKey());
	        System.out.println(" occurs " + entry.getValue() + " times");
	    }
	    
	    System.out.println("\n" + "PART 2: 10_000 latin squares of order 5 in reduced form:");
	    int[][] original5 = new int[][] { { 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 } };
	    frequencies.clear();
	    cube = createCube(original5, 5);	    
	    for ( int i = 1; i <= 10_000; i++ ) {
	        shuffleCube(cube);
	        int[][] matrix = toMatrix(cube);
	        reduce(matrix);
	        frequencies.merge(Arrays.deepToString(matrix), 1, Integer::sum);
	    }
	    
	    int count = 0;
	    for ( int frequency : frequencies.values() ) {
	        count += 1;
	        System.out.print(String.format("%s%s%2d%s%3d%s",
	        	( count > 1 ? ", " : "" ), ( count % 8 == 1 ? "\n" : "" ), count, "(", frequency, ")"));
	    }
	    
	    System.out.println("\n\n" + "PART 3: 750 latin squares of order 42, showing the last one:" + "\n");
	    int[][] matrix42 = new int[42][];
	    cube = createCube(null, 42);
	    for ( int i = 1; i <= 750; i++ ) {
	        shuffleCube(cube);
	        if ( i == 750 ) {
	        	matrix42 = toMatrix(cube);
	        	oneBased(matrix42);
	        }
	    }
	    Arrays.stream(matrix42).forEach( row -> System.out.println(Arrays.toString(row)) );
	    
	    System.out.println("\n" + "PART 4: 1,000 latin squares of order 256:" + "\n");
	    final long startTime = System.currentTimeMillis();
	    cube = createCube(null, 256);
	    for ( int i = 1; i <= 1_000; i++ ) {
	        shuffleCube(cube);
	    }
	    final long finishTime = System.currentTimeMillis();
	    System.out.println("Generated in " + ( finishTime - startTime ) + " milliseconds");
	}	
	
	private static void reduce(int[][] matrix) {
		for ( int j = 0; j < matrix.length - 1; j++ ) {
	        if ( matrix[0][j] != j ) {
	            for ( int k = j + 1; k < matrix.length; k++ ) {
	                if ( matrix[0][k] == j ) {
	                    for ( int i = 0; i < matrix.length; i++ ) {
	                    	final int temp = matrix[i][j];
	                    	matrix[i][j] = matrix[i][k];
	                    	matrix[i][k] = temp;
	                    }
	                    break;
	                }
	            }
	        }
	    }
		
		for ( int i = 1; i < matrix.length - 1; i++ ) {
	        if ( matrix[i][0] != i ) {
	            for ( int k = i + 1; k < matrix.length; k++ ) {
	                if ( matrix[k][0] == i ) {
	                    for ( int j = 0; j < matrix.length; j++ ) {
	                    	final int temp = matrix[i][j];
	                    	matrix[i][j] = matrix[k][j];
	                    	matrix[k][j] = temp;
	                    }
	                    break;
	                }
	            }
	        }
	    }		 
	}	
	
	private static int[][] toMatrix(int[][][] cube) {
	    int[][] matrix = new int[cube.length][cube.length];
	    for ( int i = 0; i < cube.length; i++ ) {
	        for ( int j = 0; j < cube.length; j++ ) {
	            for ( int k = 0; k < cube.length; k++ ) {
	                if ( cube[i][j][k] != 0 ) {
	                    matrix[i][j] = k;
	                    break;
	                }
	            }
	        }
	    }	      
	    return matrix;
	}	
	
	private static void shuffleCube(int[][][] cube) {
	    boolean proper = true;
	    
	    int rx, ry, rz;
	    do {
	    	rx = random.nextInt(0, cube.length);
	    	ry = random.nextInt(0, cube.length);
	    	rz = random.nextInt(0, cube.length);	    
	    } while ( cube[rx][ry][rz] != 0 );
	    
	    while ( true ) {
	        int ox = 0, oy = 0, oz = 0;	     

	        while ( cube[ox][ry][rz] != 1 ) {
	        	ox += 1;
	        }
	        while ( cube[rx][oy][rz] != 1 ) {
	        	oy += 1;
	        }
	        while ( cube[rx][ry][oz] != 1 ) {
	        	oz += 1;
	        }
	        
	        if ( ! proper ) {		        
		        if ( random.nextInt(2) == 0 ) {
		            ox += 1;
		            while ( cube[ox][ry][rz] != 1 ) {
			        	ox += 1;
			        }
		        }			       
		        if ( random.nextInt(2) == 0 ) {
		            oy += 1;
		            while ( cube[rx][oy][rz] != 1 ) {
			        	oy += 1;
			        }
		        }		       
		        if ( random.nextInt(2) == 0 ) {
		            oz += 1;
		            while ( cube[rx][ry][oz] != 1 ) {
			        	oz += 1;
			        }
		        }  
	        }
	        
	        cube[rx][ry][rz] += 1;
    	    cube[rx][oy][oz] += 1;
    	    cube[ox][ry][oz] += 1;
    	    cube[ox][oy][rz] += 1;

    	    cube[rx][ry][oz] -= 1;
    	    cube[rx][oy][rz] -= 1;
    	    cube[ox][ry][rz] -= 1;
    	    cube[ox][oy][oz] -= 1;
    	    
    	    if ( cube[ox][oy][oz] < 0 ) {
    	    	rx = ox; ry = oy; rz = oz;
    	        proper = false;
    	    } else { 
    	        break;
			}	            
	    }
	}
		
	private static int[][][] createCube(int[][] matrix, int size) {
	    int[][][] cube = new int[size][size][size];
	    for ( int i = 0; i < size; i++ ) {
	        cube[i] = new int[size][size];
	        for ( int j = 0; j < size; j++ ) {
	            cube[i][j] = new int[size];
	            final int k = ( matrix == null ) ? ( i + j ) % size : matrix[i][j] - 1;
	            cube[i][j][k] = 1;
	        }
	    }
	    return cube;
	}
	
	private static void oneBased(int[][] matrix) {
		for ( int i = 0; i < matrix.length; i++ ) {
	        for ( int j = 0; j < matrix.length; j++ ) {
	            matrix[i][j] += 1;
	        }
	    }
	}
		
	private static ThreadLocalRandom random = ThreadLocalRandom.current();

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

[[1, 2, 3, 4], [2, 4, 1, 3], [3, 1, 4, 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 2423 times
[[1, 2, 3, 4], [2, 1, 4, 3], [3, 4, 2, 1], [4, 3, 1, 2]] occurs 2614 times
[[1, 2, 3, 4], [2, 1, 4, 3], [3, 4, 1, 2], [4, 3, 2, 1]] occurs 2469 times

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

 1(167),  2(160),  3(180),  4(178),  5(165),  6(204),  7(187),  8(160), 
 9(203), 10(191), 11(183), 12(187), 13(162), 14(179), 15(175), 16(156), 
17(157), 18(180), 19(180), 20(174), 21(161), 22(180), 23(178), 24(181), 
25(168), 26(210), 27(175), 28(183), 29(166), 30(168), 31(183), 32(190), 
33(191), 34(189), 35(203), 36(173), 37(175), 38(179), 39(202), 40(216), 
41(189), 42(183), 43(175), 44(176), 45(155), 46(194), 47(149), 48(155), 
49(161), 50(182), 51(173), 52(180), 53(198), 54(171), 55(191), 56(169)

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

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

PART 4: 1,000 latin squares of order 256:

Generated in 681 milliseconds

Julia

Translation of: Go
const Cube = Vector{Vector{Vector{Int}}}
const Mat = Vector{Vector{Int}}

function reduced(m::Mat)
    n = length(m)
    r = deepcopy(m)
    for j in 1:n-1
        if r[1][j] != j
            for k in j+1:n
                if r[1][k] == j
                    for i in 1:n
                        r[i][j], r[i][k] = r[i][k], r[i][j]
                    end
                    break
                end
            end
        end
    end
    for i in 2:n-1
        if r[i][1] != i
            for k in i+1:n
                if r[k][1] == i
                    for j in 1:n
                        r[i][j], r[k][j] = r[k][j], r[i][j]
                    end
                    break
                end
            end
        end
    end
    return r
end

""" print matrix as small integers, no punctuation """
function print_matrix(m::Mat)
    n = length(m)
    padding = max(2, Int(ceil(log(10, n+1))) + 1)
    for i in 1:n
        for j in 1:n
            print(lpad(m[i][j], padding))
        end
        println()
    end
    println()
end

function shuffle_cube(c::Cube)
    n = length(c)
    proper = true
    rx, ry, rz = 0, 0, 0
    while true
        rx, ry, rz = rand(1:n, 3)
        c[rx][ry][rz] == 0 && break
    end
    while true
        ox = something(findfirst(i -> c[i][ry][rz] == 1, 1:n), n)
        oy = something(findfirst(i -> c[rx][i][rz] == 1, 1:n), n)
        oz = something(findfirst(i -> c[rx][ry][i] == 1, 1:n), n)
        if !proper
            rand() < 1/2 && (ox = something(findlast(i -> c[i][ry][rz] == 1, 1:n), n))
            rand() < 1/2 && (oy = something(findlast(i -> c[rx][i][rz] == 1, 1:n), n))
            rand() < 1/2 && (oz = something(findlast(i -> c[rx][ry][i] == 1, 1:n), n))
        end

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

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

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

function matrix(c::Cube)::Mat
    n = length(c)
    m = [[0 for i in 1:n] for j in 1:n]
    for i in 1:n, j in 1:n
        for k in 1:n
            if c[i][j][k] != 0
                m[i][j] = k
                break
            end
        end
    end
    return m
end

function cube(from, n)
    c = [[[0 for i in 1:n] for j in 1:n] for k in 1:n]
    for i in 1:n, j in 1:n
        k = (from isa Nothing) ? mod1(i + j, n) : from[i][j]
        c[i][j][k] = 1
    end
    return c
end

function testJacobsenMatthews()
    # part 1
    println("PART 1: 10,000 latin Squares of order 4 in reduced form:\n")
    from = [[1, 2, 3, 4], [2, 1, 4, 3], [3, 4, 1, 2], [4, 3, 2, 1]]
    freqs4 = Dict{Array, Int}()
    c = cube(from, 4)
    for i in 1:10000
        shuffle_cube(c)
        m = matrix(c)
        rm = reduced(m)
        n = get!(freqs4, rm, 0)
        freqs4[rm] = n + 1
    end
    for (a, freq) in freqs4
        print_matrix(a)
        println("Occurs $freq times\n")
    end

    # part 2
    println("\nPART 2: 10,000 latin squares of order 5 in reduced form:\n")
    from = [[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 = Dict{Array, Int}()
    c = cube(from, 5)
    for i in 1:10000
        shuffle_cube(c)
        m = matrix(c)
        rm = reduced(m)
        n = get!(freqs5, rm, 0)
        freqs5[rm] = n + 1
    end
    for (i, freq) in enumerate(sort(collect(values(freqs5))))
        i > 1 && (print(", "); (i - 1) % 8 == 0 && println())
        print(lpad(i, 2), "(", lpad(freq, 3), ")")
    end
    println("\n")

    # part 3
    println("\nPART 3: 750 latin squares of order 42, showing the last one:\n")
    m42 = [[0 for i in 1:42] for j in 1:42]
    c = cube(nothing, 42)
    for i in 1:750
        shuffle_cube(c)
        i == 750 && (m42 = matrix(c))
    end
    print_matrix(m42)

    # part 4
    println("\nPART 4: 1000 latin squares of order 256:\n")
    @time begin
        c = cube(nothing, 256)
        for i in 1:1000
            shuffle_cube(c)
        end
    end
end

testJacobsenMatthews()
Output:
PART 1: 10,000 latin Squares of order 4 in reduced form:

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

Occurs 2508 times

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

Occurs 2427 times

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

Occurs 2529 times

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

Occurs 2536 times


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

 1(152),  2(152),  3(153),  4(154),  5(158),  6(160),  7(160),  8(160), 
 9(162), 10(165), 11(166), 12(167), 13(168), 14(170), 15(170), 16(172),
17(172), 18(173), 19(174), 20(174), 21(175), 22(177), 23(177), 24(177),
25(177), 26(178), 27(179), 28(180), 29(180), 30(181), 31(181), 32(182),
33(182), 34(182), 35(183), 36(184), 37(185), 38(185), 39(185), 40(186), 
41(187), 42(187), 43(187), 44(188), 45(189), 46(189), 47(190), 48(195),
49(195), 50(197), 51(197), 52(199), 53(199), 54(199), 55(201), 56(203)


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

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


PART 4: 1000 latin squares of order 256:

 10.811605 seconds (745.43 k allocations: 157.305 MiB, 0.30% gc time)

Nim

Translation of: Go
import random, sequtils, strformat, tables, times

type
  Vector = seq[int]
  Matrix = seq[Vector]
  Cube = seq[Matrix]
  Array16 = array[16, int]
  Array25 = array[25, int]


func newCube(m: Matrix; n: Positive): Cube =
  result.setLen(n)
  for i in 0..<n:
    result[i].setLen(n)
    for j in 0..<n:
      result[i][j].setLen(n)
      let k = if m.len == 0: (i + j) mod n else: m[i][j] - 1
      result[i][j][k] = 1


proc shuffle(c: var Cube) =
  let n = c[0].len
  var proper = true

  var rx, ry, rz: int
  while true:
    rx = rand(n - 1)
    ry = rand(n - 1)
    rz = rand(n - 1)
    if c[rx][ry][rz] == 0: break

  while true:
    var ox, oy, oz = 0

    while ox < n:
      if c[ox][ry][rz] == 1: break
      inc ox
    if not proper and rand(1) == 0:
      inc ox
      while ox < n:
        if c[ox][ry][rz] == 1: break
        inc ox

    while oy < n:
      if c[rx][oy][rz] == 1: break
      inc oy
    if not proper and rand(1) == 0:
      inc oy
      while oy < n:
        if c[rx][oy][rz] == 1: break
        inc oy

    while oz < n:
      if c[rx][ry][oz] == 1: break
      inc oz
    if not proper and rand(1) == 0:
      inc oz
      while oz < n:
        if c[rx][ry][oz] == 1: break
        inc oz

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

    dec c[rx][ry][oz]
    dec c[rx][oy][rz]
    dec c[ox][ry][rz]
    dec 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 =
  let n = c[0].len
  result = newSeqWith(n, newSeq[int](n))
  for i in 0..<n:
    for j in 0..<n:
      for k in 0..<n:
        if c[i][j][k] != 0:
          result[i][j] = k
          break


func toReduced(m: Matrix): Matrix =
  let n = m.len
  result = m

  for j in 0..n-2:
    if result[0][j] != j:
      for k in j+1..<n:
        if result[0][k] == j:
          for i in 0..<n:
            swap result[i][j], result[i][k]
          break

  for i in 1..n-2:
    if result[i][0] != i:
      for k in i+1..<n:
        if result[k][0] == i:
          for j in 0..<n:
            swap result[i][j], result[k][j]
          break


func asArray16(m: Matrix): Array16 =
  var k = 0
  for i in 0..3:
    for j in 0..3:
      result[k] = m[i][j]
      inc k


func asArray25(m: Matrix): Array25 =
  var k = 0
  for i in 0..4:
    for j in 0..4:
      result[k] = m[i][j]
      inc k


proc printArray16(a: Array16) =
  for i in 0..3:
    for j in 0..3:
      let k = i * 4 + j
      stdout.write &"{a[k]+1:2} "   # Back to 1 based.
    echo()
  echo()


proc printMatrix(m: Matrix) =
  let n = m.len
  for i in 0..<n:
    for j in 0..<n:
      stdout.write &"{m[i][j]+1:2} "  # Back to 1 based.
    echo()
  echo()


randomize()

# Part 1.
echo "Part 1: 10_000 latin Squares of order 4 in reduced form:\n"
const From1: Matrix = @[@[1, 2, 3, 4], @[2, 1, 4, 3], @[3, 4, 1, 2], @[4, 3, 2, 1]]
var freqs4: CountTable[Array16]
var c = newCube(From1, 4)
for _ in 1..10_000:
  c.shuffle()
  let m = c.toMatrix
  let rm = m.toReduced
  let a16 = rm.asArray16
  freqs4.inc(a16)

for a, freq in freqs4.pairs:
  printArray16(a)
  echo &"Occurs {freq} times\n"

# Part 2.
echo "\nPart 2: 10_000 latin squares of order 5 in reduced form:"
const From2: 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]]
var freqs5: CountTable[Array25]
c = newCube(From2, 5)
for _ in 1..10_000:
  c.shuffle()
  let m = c.toMatrix
  let rm = m.toReduced
  let a25 = rm.asArray25
  freqs5.inc(a25)

var count = 0
for freq in freqs5.values:
  inc count
  if count > 1: stdout.write ", "
  if (count - 1) mod 8 == 0: echo()
  stdout.write &"{count:2}({freq:3})"
echo '\n'

# Part 3.
echo "\nPart 3: 750 latin squares of order 42, showing the last one:\n"
var m42: Matrix
c = newCube(@[], 42)
for i in 1..750:
  c.shuffle()
  if i == 750:
    m42 = c.toMatrix
printMatrix(m42)

# Part 4.
echo "\nPart 4: 1000 latin squares of order 256:\n"
let start = cpuTime()
c = newCube(@[], 256)
for _ in 1..1000:
  c.shuffle()
echo &"Generated in {cpuTime() - start:.3f} s."
Output:
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 2469 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  3  4  1 
 3  4  1  2 
 4  1  2  3 

Occurs 2561 times

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

Occurs 2540 times


Part 2: 10_000 latin squares of order 5 in reduced form:

 1(175),  2(201),  3(209),  4(158),  5(175),  6(190),  7(178),  8(173), 
 9(172), 10(179), 11(199), 12(170), 13(167), 14(199), 15(166), 16(164), 
17(180), 18(189), 19(170), 20(183), 21(165), 22(172), 23(199), 24(190), 
25(165), 26(178), 27(178), 28(190), 29(173), 30(187), 31(163), 32(157), 
33(179), 34(190), 35(168), 36(209), 37(185), 38(181), 39(161), 40(180), 
41(162), 42(161), 43(184), 44(157), 45(195), 46(179), 47(181), 48(161), 
49(198), 50(190), 51(176), 52(184), 53(148), 54(209), 55(167), 56(181)


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

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


Part 4: 1000 latin squares of order 256:

Generated in 1.306 s.

Phix

Translation of: Go
function shuffleCube(sequence c)
    integer n = length(c), rx, ry, rz
    bool bProper = true
    while true do
        rx = rand(n)
        ry = rand(n)
        rz = rand(n)
        if c[rx][ry][rz] == 0 then exit end if
    end while
    while true do
        integer ox, oy, oz
        for i=1 to n do
            ox = i
            if c[ox][ry][rz] == 1 then exit end if
        end for
        if not bProper and rand(2)==2 then
            for i=ox+1 to n do
                ox = i
                if c[ox][ry][rz] == 1 then exit end if
            end for
        end if
        for i=1 to n do
            oy = i
            if c[rx][oy][rz] == 1 then exit end if
        end for
        if not bProper and rand(2)==2 then
            for i=oy+1 to n do
                oy = i
                if c[rx][oy][rz] == 1 then exit end if 
            end for
        end if
        for i=1 to n do
            oz = i
            if c[rx][ry][oz] == 1 then exit end if
        end for
        if not bProper and rand(2)==2 then
            for i=oz+1 to n do
                oz = i
                if c[rx][ry][oz] == 1 then exit end if
            end for
        end if
 
        c[rx][ry][rz] += 1
        c[rx][oy][oz] += 1
        c[ox][ry][oz] += 1
        c[ox][oy][rz] += 1
 
        c[rx][ry][oz] -= 1
        c[rx][oy][rz] -= 1
        c[ox][ry][rz] -= 1
        c[ox][oy][oz] -= 1
 
        if c[ox][oy][oz] < 0 then
            {rx, ry, rz} = {ox, oy, oz}
            bProper = false
        else
            bProper = true
            exit
        end if
    end while
    return c
end function
 
function toMatrix(sequence c)
    integer n = length(c)
    sequence m = repeat(repeat(0,n),n)
    for i=1 to n do
        for j=1 to n do
            for k=1 to n do
                if c[i][j][k] != 0 then
                    m[i][j] = k
                    exit
                end if
            end for
        end for
    end for
    return m
end function
 
function toReduced(sequence m)
    integer n := length(m)
    m = deep_copy(m)
    for j=1 to n-1 do
        if m[1][j]!=j then
            for k=j+1 to n do
                if m[1][k]==j then
                    for i=1 to n do
                        atom mij = m[i][j],
                             mik = m[i][k]
                        m[i][j] = mik
                        m[i][k] = mij

                    end for
                    exit
                end if
            end for
        end if
    end for
    for i=2 to n-1 do
        if m[i][1]!=i then
            for k=i+1 to n do
                if m[k][1]==i then
                    for j=1 to n do
                        atom mij = m[i][j],
                             mkj = m[k][j]
                        m[i][j] = mkj
                        m[k][j] = mij
                    end for
                    exit
                end if
            end for
        end if
    end for
    return m
end function
 
function makeCube(object orig, integer n)
    sequence c = repeat(repeat(repeat(0,n),n),n)
    for i=1 to n do
        for j=1 to n do
            integer k = iff(orig==NULL?mod(i+j,n)+1:orig[i][j])
            c[i][j][k] = 1
        end for
    end for
    return c
end function
 
procedure main()
 
    printf(1,"Part 1: 10,000 latin Squares of order 4 in reduced form:\n\n")
    sequence orig = {{1, 2, 3, 4}, {2, 1, 4, 3}, {3, 4, 1, 2}, {4, 3, 2, 1}},
             c := makeCube(orig, 4), m, rm, fk
    integer freq = new_dict()
    for i=1 to 10000 do
        c = shuffleCube(c)
        m = toMatrix(c)
        rm = toReduced(m)
        setd(rm,getd(rm,freq)+1,freq)
    end for
    fk = getd_all_keys(freq)
    for i=1 to length(fk) do
        printf(1,"%v occurs %d times\n", {fk[i],getd(fk[i],freq)})
    end for
 
    printf(1,"\nPart 2: 10,000 latin squares of order 5 in reduced form:\n\n")
    orig = {{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}}
    c = makeCube(orig, 5)
    destroy_dict(freq, justclear:=true)
    for i=1 to 10000 do
        c = shuffleCube(c)
        m = toMatrix(c)
        rm = toReduced(m)
        setd(rm,getd(rm,freq)+1,freq)
    end for
    fk = getd_all_keys(freq)
    for i=1 to length(fk) do
        fk[i] = sprintf("%2d(%3d)", {i,getd(fk[i],freq)})
    end for
    puts(1,join_by(fk,8,7," ","\n"))
    destroy_dict(freq)
 
    -- part 3
    printf(1,"\nPart 3: 750 latin squares of order 42, showing the last one:\n\n")
    c = makeCube(NULL, 42)
    for i=1 to 750 do
        c = shuffleCube(c)
    end for
    m = toMatrix(c)
    integer n := length(m)
    for i=1 to n do
        for j=1 to n do
            m[i,j] = sprintf("%2d",m[i,j])
        end for
        m[i] = join(m[i]," ")
    end for
    printf(1,"%s\n",join(m,"\n"))
 
    -- part 4
    printf(1,"\nPART 4: 1000 latin squares of order 256:\n\n")
    atom t0 = time()
    c = makeCube(NULL, 256)
    for i=1 to 1000 do
        c = shuffleCube(c)
    end for
    printf(1,"Generated in %s\n", elapsed(time()-t0))
end procedure
main()
Output:
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 2503 times
{{1,2,3,4},{2,1,4,3},{3,4,2,1},{4,3,1,2}} occurs 2560 times
{{1,2,3,4},{2,3,4,1},{3,4,1,2},{4,1,2,3}} occurs 2510 times
{{1,2,3,4},{2,4,1,3},{3,1,4,2},{4,3,2,1}} occurs 2427 times

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

 1(172)  9(197) 17(228) 25(166) 33(171) 41(224) 49(171)
 2(168) 10(162) 18(216) 26(227) 34(172) 42(155) 50(226)
 3(159) 11(198) 19(206) 27(165) 35(189) 43(190) 51(174)
 4(170) 12(207) 20(159) 28(166) 36(177) 44(171) 52(196)
 5(211) 13(148) 21(172) 29(173) 37(183) 45(189) 53(197)
 6(169) 14(163) 22(128) 30(179) 38(184) 46(138) 54(173)
 7(168) 15(155) 23(146) 31(170) 39(187) 47(170) 55(206)
 8(193) 16(177) 24(146) 32(176) 40(157) 48(183) 56(177)

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

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

PART 4: 1000 latin squares of order 256:

Generated in 19.5s

Unfortunately the last part of this task exposes the relatively poor performance of subscripting in phix.

Raku

Translation of: Go
# 20210729 Raku programming solution 

#!/usr/bin/env raku

sub makeCube(\from, Int \n) {  
   my @c = [[[ 0 xx n ] xx n ] xx n ]; 
   from.Bool ?? do race for ^n X ^n -> (\i,\j) { @c[i;j; { from[i;j]-1 } ] = 1 }
             !! do race for ^n X ^n -> (\i,\j) { @c[i;j; {   (i+j)%n   } ] = 1 }
   return @c                   
}

sub shuffleCube(@c) {
   my ($rx, $ry, $rz); my \n = +@c; my Bool \proper = $ = True;

   repeat { ($rx ,$ry, $rz) = (^n).roll: 3 } until @c[$rx;$ry;$rz] == 0;
   loop { 
      my ($ox, $oy, $oz); 
      for ^n { last if @c[ $ox = $_ ;$ry;$rz] == 1 }
      if !proper and (^3).roll==0 { 
         for $ox^…^n { last if @c[ $ox = $_ ;$ry;$rz] == 1 } 
      }
      for ^n { last if @c[$rx; $oy = $_ ;$rz] == 1 }
      if !proper and (^3).roll==0 {
         for $oy^…^n { last if @c[$rx; $oy = $_ ;$rz] == 1 }
      } 
      for ^n { last if @c[$rx;$ry;  $oz = $_ ] == 1 }
      if !proper and (^3).roll==0 {
         for $oz^…^n { last if @c[$rx;$ry; $oz = $_ ] == 1 }
      }

      (@c[$rx;$ry;$rz],@c[$rx;$oy;$oz],@c[$ox;$ry;$oz],@c[$ox;$oy;$rz]) »+=»1;
      (@c[$rx;$ry;$oz],@c[$rx;$oy;$rz],@c[$ox;$ry;$rz],@c[$ox;$oy;$oz]) »-=»1; 
      
      @c[$ox;$oy;$oz] < 0 ?? (($rx,$ry,$rz) = ($ox,$oy,$oz)) !! last ;
      proper = False 
   }        
}

sub toMatrix(@c) {
   my \n = +@c;  
   my @m = [[0 xx n] xx n]; 
   for ^n X ^n -> (\i,\j) { 
      for ^n -> \k { if @c[i;j;k] != 0 { @m[i;j] = k and last } }
   }
   return @m
}

sub toReduced(@m is copy) { 
   my \n = +@m; 
   for 0…(n-2) -> \j {
      if ( @m[0;j] != j ) { 
         for j^…^n -> \k {
            if ( @m[0;k] == j ) { 
               for 0…^n -> \i { (@m[i;j], @m[i;k]) = (@m[i;k], @m[i;j]) }
               last 
            } 
         } 
      } 
   }   
   for 1…(n-2) -> \i {
      if ( @m[i;0] != i ) {
         for i^…^n -> \k { 
            if ( @m[k;0] == i ) { 
               for 0…^n -> \j { (@m[i;j], @m[k;j]) = (@m[k;j], @m[i;j]) }     
               last      
            }
         }
      } 
   } 
   return @m
}

sub printAs1based { say ($_ »+» 1).Str for @_.rotor: @_.elems.sqrt }

my (%freq, @c, @in);

say "Part 1: 10,000 latin Squares of order 4 in reduced form:\n";
@in = [[1, 2, 3, 4], [2, 1, 4, 3], [3, 4, 1, 2], [4, 3, 2, 1]];
@c = makeCube(@in, 4); 
for ^10_000 {
   shuffleCube @c ; 
   %freq{@c.&toMatrix.&toReduced».List.flat.Str}++
}    
for %freq.kv -> $k, $v { 
   printAs1based $k.split(' ');   
   say "\nOccurs $v times.\n"
}

say "Part 2: 10,000 latin Squares of order 5 in reduced form:\n";
@in = [ [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] ];
%freq = ();
@c = makeCube(@in, 5);
for ^10_000 { 
   shuffleCube @c ; 
   %freq{@c.&toMatrix.&toReduced».List.flat.Str}++ 
}
for %freq.values.kv -> $i, $j { printf "%2d(%3d)%s", $i+1, $j, ' ' }

say "\n\nPart 3: 750 latin squares of order 42, showing the last one:\n";
@c = makeCube([], 42); # (1..42).pick(*)
( for ^750 { shuffleCube @c } ) and printAs1based @c.&toMatrix».List.flat ;

say "\nPart 4: 100 latin squares of order 256:\n";
my $snapshot = now;
@c = makeCube([], 256);
for ^100 { shuffleCube @c } # without hyper, will do only 100 cycles 
say "Generated in { now - $snapshot } seconds."
Output:
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 2442 times.

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

Occurs 2705 times.

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

Occurs 2548 times.

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

Occurs 2305 times.

Part 2: 10,000 latin Squares of order 5 in reduced form:

 1(210)  2(191)  3(186)  4(158)  5(219)  6(164)  7(147)  8(160)  9(196) 10(188) 11(193) 12(168) 13(195) 14(173) 15(151) 16(184) 17(211) 18(171) 19(185) 20(155) 21(157) 22(191) 23(195) 24(177) 25(157) 26(191) 27(165) 28(178) 29(191) 30(180) 31(193) 32(176) 33(196) 34(178) 35(156) 36(168) 37(155) 38(152) 39(155) 40(223) 41(159) 42(165) 43(210) 44(175) 45(195) 46(188) 47(178) 48(154) 49(172) 50(176) 51(178) 52(164) 53(175) 54(196) 55(187) 56(189) 

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

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

Part 4: 100 latin squares of order 256:

Generated in 76.816295878 seconds.

Wren

Translation of: Go
Library: Wren-fmt
Library: Wren-seq
import "random" for Random
import "./fmt" for Fmt
import "./seq" for Lst

var rand = Random.new()

var toReduced = Fn.new { |m|
    var n = m.count
    var r = List.filled(n, null)
    for (i in 0...n) r[i] = m[i].toList
    for (j in 0...n-1) {
        if (r[0][j] != j) {
            for (k in j+1...n) {
                if (r[0][k] == j) {
                    for (i in 0...n) r[i].swap(j, k)
                    break
                }
            }
        }
    }
    for (i in 1...n-1) {
        if (r[i][0] != i) {
            for (k in i+1...n) {
                if (r[k][0] == i) {
                    for (j in 0...n) {
                        var t = r[i][j]
                        r[i][j] = r[k][j]
                        r[k][j] = t
                    }
                    break
                }
            }
        }
    }
    return r
}

// 'm' is assumed to be 0 based
var printMatrix = Fn.new { |m|
    var n = m.count
    for (i in 0...n) {
        for (j in 0...n) Fmt.write("$2d ", m[i][j]+1) // back to 1 based
        System.print()
    }
    System.print()
}

// converts 4 x 4 matrix to 'flat' list
var asList16 = Fn.new { |m| Lst.flatten(m) }

// converts 5 x 5 matrix to 'flat' list
var asList25 = Fn.new { |m| Lst.flatten(m) }

// 'a' is assumed to be 0 based
var printList16 = Fn.new { |a|
    for (i in 0...4) {
        for (j in 0...4) {
            var k = i*4 + j
            Fmt.write("$2d ", a[k]+1) // back to 1 based
        }
        System.print()
    }
    System.print()
}

var shuffleCube = Fn.new { |c|
    var n = c[0].count
    var proper = true
    var rx
    var ry
    var rz
    while (true) {
        rx = rand.int(n)
        ry = rand.int(n)
        rz = rand.int(n)
        if (c[rx][ry][rz] == 0) break
    }
    while (true) {
        var ox = 0
        var oy = 0
        var oz = 0
        while (ox < n) {
            if (c[ox][ry][rz] == 1) break
            ox = ox + 1
        }
        if (!proper && rand.int(2) == 0) {
            ox = ox + 1
            while (ox < n) {
                if (c[ox][ry][rz] == 1) break
                ox = ox + 1
            }
        }
        while (oy < n) {
            if (c[rx][oy][rz] == 1) break
            oy = oy + 1
        }
        if (!proper && rand.int(2) == 0) {
            oy = oy + 1
            while (oy < n) {
                if (c[rx][oy][rz] == 1) break
                oy = oy + 1
            }
        }
        while (oz < n) {
            if (c[rx][ry][oz] == 1) break
            oz = oz + 1
        }
        if (!proper && rand.int(2) == 0) {
            oz = oz + 1
            while (oz < n) {
                if (c[rx][ry][oz] == 1) break
                oz = oz + 1
            }
        }

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

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

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

var toMatrix = Fn.new { |c|
    var n = c[0].count
    var m = List.filled(n, null)
    for (i in 0...n) m[i] = List.filled(n, 0)
    for (i in 0...n) {
        for (j in 0...n) {
            for (k in 0...n) {
                if (c[i][j][k] != 0) {
                    m[i][j] = k
                    break
                }
            }
        }
    }
    return m
}

// 'from' matrix is assumed to be 1 based
var makeCube = Fn.new { |from, n|
    var c = List.filled(n, null)
    for (i in 0...n) {
        c[i] = List.filled(n, null)
        for (j in 0...n) {
            c[i][j] = List.filled(n, 0)
            var k = (!from) ? (i + j) % n : from[i][j] - 1
            c[i][j][k] = 1
        }
    }
    return c
}

// part 1
System.print("PART 1: 10,000 latin Squares of order 4 in reduced form:\n")
var from = [ [1, 2, 3, 4], [2, 1, 4, 3], [3, 4, 1, 2], [4, 3, 2, 1] ]
var freqs4 = {}
var c = makeCube.call(from, 4)
for (i in 1..10000) {
    shuffleCube.call(c)
    var m = toMatrix.call(c)
    var rm = toReduced.call(m)
    var a16 = asList16.call(rm)
    var a16s = a16.toString // can't use a list as a map key so convert it to string
    freqs4[a16s] = freqs4[a16s] ? freqs4[a16s] + 1 : 1
}
for (me in freqs4) {
    printList16.call(me.key[1..-2].split(", ").map { |n| Num.fromString(n) }.toList)
    Fmt.print("Occurs $d times\n", me.value)
}

// part 2
System.print("\nPART 2: 10,000 latin squares of order 5 in reduced form:")
from = [ [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] ]
var freqs5 = {}
c = makeCube.call(from, 5)
for (i in 1..10000) {
    shuffleCube.call(c)
    var m = toMatrix.call(c)
    var rm = toReduced.call(m)
    var a25 = asList25.call(rm)
    var a25s = a25.toString // can't use a list as a map key so convert it to string
    freqs5[a25s] = freqs5[a25s] ? freqs5[a25s] + 1 : 1
}
var count = 0
for (freq in freqs5.values) {
    count = count + 1
    if (count > 1) System.write(", ")
    if ((count-1) % 8 == 0) System.print()
    Fmt.write("$2d($3d)", count, freq)
}
System.print("\n")

// part 3
System.print("\nPART 3: 750 latin squares of order 42, showing the last one:\n")
var m42
c = makeCube.call(null, 42)
for (i in 1..750) {
    shuffleCube.call(c)
    if (i == 750) m42 = toMatrix.call(c)
}
printMatrix.call(m42)

// part 4
System.print("\nPART 4: 1,000 latin squares of order 256:\n")
var start = System.clock
c = makeCube.call(null, 256)
for (i in 1..1000) shuffleCube.call(c)
var elapsed = System.clock - start
Fmt.print("Generated in $s seconds", elapsed)
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 2510 times

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

Occurs 2498 times

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

Occurs 2506 times

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

Occurs 2486 times


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

 1(187),  2(179),  3(186),  4(176),  5(177),  6(189),  7(193),  8(182), 
 9(168), 10(169), 11(147), 12(200), 13(198), 14(169), 15(200), 16(173), 
17(184), 18(179), 19(151), 20(174), 21(160), 22(198), 23(153), 24(184), 
25(170), 26(180), 27(171), 28(180), 29(184), 30(178), 31(197), 32(173), 
33(185), 34(181), 35(200), 36(188), 37(176), 38(196), 39(193), 40(183), 
41(163), 42(163), 43(173), 44(178), 45(177), 46(160), 47(155), 48(165), 
49(181), 50(188), 51(187), 52(182), 53(162), 54(192), 55(183), 56(180)


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

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


PART 4: 1,000 latin squares of order 256:

Generated in 10.828862 seconds
Cookies help us deliver our services. By using our services, you agree to our use of cookies.