1## Copyright (C) 2009 Jaroslav Hajek <highegg@gmail.com> 2## 3## This program is free software; you can redistribute it and/or modify it under 4## the terms of the GNU General Public License as published by the Free Software 5## Foundation; either version 3 of the License, or (at your option) any later 6## version. 7## 8## This program is distributed in the hope that it will be useful, but WITHOUT 9## ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or 10## FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more 11## details. 12## 13## You should have received a copy of the GNU General Public License along with 14## this program; if not, see <http://www.gnu.org/licenses/>. 15 16## -*- texinfo -*- 17## @deftypefn{Function File} [@var{x}, @var{ntrial}] = solvesudoku (@var{s}) 18## Solves a classical 9x9 sudoku. @var{s} should be a 9x9 array with 19## numbers from 0:9. 0 indicates empty field. 20## Returns the filled table or empty matrix if no solution exists. 21## If requested, @var{ntrial} returns the number of trial-and-error steps needed. 22## @end deftypefn 23 24## This uses a recursive backtracking technique combined with revealing new singleton 25## fields by logic. The beauty of it is that it is completely vectorized. 26 27function [x, ntrial] = solvesudoku (s) 28 29 if (nargin != 1) 30 print_usage (); 31 endif 32 33 if (! (ismatrix (s) && ndims (s) == 2 && all (size (s) == [9, 9]))) 34 error ("needs a 9x9 matrix"); 35 endif 36 37 if (! ismember (unique (s(:)), 0:9)) 38 error ("matrix must contain values from 0:9"); 39 endif 40 41 if (! verifysudoku (s)) 42 error ("matrix is not a valid sudoku grid"); 43 endif 44 45 [x, ntrial] = solvesudoku_rec (s); 46 47endfunction 48 49function ok = verifysudoku (s) 50 [i, j, k] = find (s); 51 b = false (9, 9, 9); 52 b(sub2ind ([9, 9, 9], i, j, k)) = true; 53 okc = sum (b, 1) <= 1; 54 okr = sum (b, 2) <= 1; 55 b = reshape (b, [3, 3, 3, 3, 9]); 56 ok3 = sum (sum (b, 1), 3) <= 1; 57 ok = all (okc(:) & okr(:) & ok3(:)); 58endfunction 59 60function [x, ntrial] = solvesudoku_rec (s) 61 62 x = s; 63 ntrial = 0; 64 65 ## Run until the logic is exhausted. 66 do 67 b = getoptions (x); 68 s = x; 69 x = getsingletons (b, x); 70 finished = isempty (x) || all (x(:)); 71 until (finished || all ((x == s)(:))); 72 73 if (! finished) 74 x = []; 75 ## Find the field with minimum possibilities. 76 sb = sum (b, 3); 77 sb(s != 0) = 10; 78 [msb, i] = min (sb(:)); 79 [i, j] = ind2sub ([9, 9], i); 80 ## Try all guesses. 81 for k = find (b(i,j,:))' 82 s(i,j) = k; 83 [x, ntrial1] = solvesudoku_rec (s); 84 ntrial += 1 + ntrial1; 85 if (! isempty (x)) 86 ## Found solutions. 87 break; 88 endif 89 s(i,j) = 0; 90 endfor 91 endif 92 93endfunction 94 95## Given a 9x9x9 logical array of allowed values, get the logical singletons. 96function s = getsingletons (b, s) 97 98 n0 = sum (s(:) != 0); 99 100 ## Check for fields with only one option. 101 sb = sum (b, 3); 102 if (any (sb(:) == 0)) 103 s = []; 104 return; 105 else 106 s1 = sb == 1; 107 ## We want to return as soon as some new singletons are found. 108 [s(s1), xx] = find (reshape (b, [], 9)(s1, :).'); 109 if (sum (s(:) != 0) > n0) 110 return; 111 endif 112 endif 113 114 ## Check for columns where a number has only one field left. 115 sb = squeeze (sum (b, 1)); 116 if (any (sb(:) == 0)) 117 s = []; 118 return; 119 else 120 s1 = sb == 1; 121 [j, k] = find (s1); 122 [i, xx] = find (b(:, s1)); 123 s(sub2ind ([9, 9], i, j)) = k; 124 if (sum (s(:) != 0) > n0) 125 return; 126 endif 127 endif 128 129 ## Ditto for rows. 130 sb = squeeze (sum (b, 2)); 131 if (any (sb(:) == 0)) 132 s = []; 133 return; 134 else 135 s1 = sb == 1; 136 [i, k] = find (s1); 137 [j, xx] = find (permute (b, [2, 1, 3])(:, s1)); 138 s(sub2ind ([9, 9], i, j)) = k; 139 if (sum (s(:) != 0) > n0) 140 return; 141 endif 142 endif 143 144 ## 3x3 tiles. 145 bb = reshape (b, [3, 3, 3, 3, 9]); 146 sb = squeeze (sum (sum (bb, 1), 3)); 147 if (any (sb(:) == 0)) 148 s = []; 149 return; 150 else 151 s1 = reshape (sb == 1, 9, 9); 152 [j, k] = find (s1); 153 [i, xx] = find (reshape (permute (bb, [1, 3, 2, 4, 5]), 9, 9*9)(:, s1)); 154 [i1, i2] = ind2sub ([3, 3], i); 155 [j1, j2] = ind2sub ([3, 3], j); 156 s(sub2ind ([3, 3, 3, 3], i1, j1, i2, j2)) = k; 157 if (sum (s(:) != 0) > n0) 158 return; 159 endif 160 endif 161 162endfunction 163 164## Given known values (singletons), calculate options. 165function b = getoptions (s) 166 167 ## Find true values. 168 [i, j, s] = find (s); 169 ## Columns. 170 bc = true (9, 9, 9); 171 bc(:, sub2ind ([9, 9], j, s)) = false; 172 ## Rows. 173 br = true (9, 9, 9); 174 br(:, sub2ind ([9, 9], i, s)) = false; 175 ## 3x3 tiles. 176 b3 = true (3, 3, 3, 3, 9); 177 b3(:, :, sub2ind ([3, 3, 9], ceil (i/3), ceil (j/3), s)) = false; 178 ## Permute elements to correct order. 179 br = permute (br, [2, 1, 3]); 180 b3 = reshape (permute (b3, [1, 3, 2, 4, 5]), [9, 9, 9]); 181 ## The singleton fields themselves. 182 bb = true (9*9, 9); 183 bb(sub2ind ([9, 9], i, j), :) = false; 184 bb = reshape (bb, [9, 9, 9]); 185 ## Form result. 186 b = bc & br & b3 & bb; 187 ## Correct singleton fields. 188 b = reshape (b, 9, 9, 9); 189 b(sub2ind ([9, 9, 9], i, j, s)) = true; 190 191endfunction 192 193