探究一個問題:全體數獨題的總解數是多少?
問題:全體數獨題的總解數是多少?
這個問題,在最初寫 SudokuSolver 1.0 時,就已經想到了。在 SudokuSolver 1.0:用C++實現的數獨解題程式 【二】 的末尾專門提到了全空數獨題,並指出:它的解是最多的,事實上任意一個數獨題的解都是它的解。
這裡要探究的問題實際就是全空數獨題有多少個解的問題。前兩天在 SudokuSolver 2.4 的基礎上實現了一個 runrrun <sum> 的命令,比如 runrun 10000 是讓程式接到命令開始找到一萬個新解(或把找到全體剩餘的解)才停止回到交互狀態準備接受下一個命令。這個新版的程式陸續跑了幾天,最新取樣的結果是:
Order please: runrun 41000000 ... Run-run time: 5749367 milliseconds; current solutions: 41000000 steps: 0 # 1 # 120797048 total solutions: 0 # 700123572. Order please: show 123 456 789 456 789 123 789 123 456 245 867 931 967 312 548 831 945 672 572 698 314 618 534 297 394 271 865 Fine [steps:120797048, solution sum:700123572]. Order please:
從輸出的資訊看,當時求得的總解數已經超過 7 億,而且這些解的前三行填值完全一樣,即:
123 456 789 456 789 123 789 123 456
由此也可大致推斷最終的總解數會是一個非常大的數,為說明方便,把這個數記作 δ。
另外,還可以看到 steps 統計值發生了溢出,即超出了 232 被截斷為 120797048。這在意料之中,步數比解數增長更快,更容易溢出。為了應對溢出問題,改進的程式在 m_steps 基礎上又增加了兩個分量:m_bigSteps 和 m_bigbigSteps,因此總步數會等於 m_bigbigSteps × 264 + m_bigSteps × 232 + m_steps;同樣,在 m_soluSum 基礎上增加了 m_bigSum,總解數會等於 m_bigSum × 232 + m_soluSum。
δ 的上界估算
在進一步探討改進程式之前,先用數學方法大致估算一下 δ 的上界。
假設全體解中第一行填值為 123 456 789 的總解數為 S。記 A = {1,2,3,4,5,6,7,8,9},則 A 到自身的雙射會有 9! 種,考慮其中任意一個非恆等變換,比如 (1,2,3,4,5,6,7,8,9) -> (2,3,4,5,6,7,8,1),會得到 S 個第一行填值為 234 567 891 的解,因此 δ = 9! · S。
9! = 362,880,假設用程式對如下的數獨題(後面簡記為 1R 數獨題)
123 456 789 000 000 000 000 000 000 000 000 000 000 000 000 000 000 000 000 000 000 000 000 000 000 000 000
求其總解數 S 需要一年的時間,那麼對全空數獨題用程式求解其總解數 δ 則需要 9! = 362,880 年。
繼續估算 S 的上界。在 1R 數獨題中,逐個考察第 1 宮中的 6 個空位:[2,1] 位置上有 6 個候選值;選定 [2,1] 位置取值後,[2,2] 位置上有 5 個候選值;依此類推,即如下所示:
123 456 789 654 000 000 321 000 000
這裡加粗的數字不是表示對應空位上的填值,而是對應空位上填值的可能個數。這說明第 1 宮裡的 6 個空位的整體填值會有 6! 種可能。同樣,考察第 1 列,又會得到如下的情形:
123 456 789 654 000 000 321 000 000 600 000 000 500 000 000 ① 400 000 000 300 000 000 200 000 000 100 000 000
即第 1 列的第 2 節和第 3 節構成的 6 個空位整體會有 6! 種填值可能。
繼續依次考察第 2 行的 6 個空位,[2,4] 位置上的值不能為 4、5、6,也不能為第 2 行第一節里的三個數 a、b、c,但這時會出現不確定的情形,比如 a、b、c 都選自 {4,5,6},則 [2,4] 位置上的值會有 6 種可能,而如果 a、b、c 都選自 {7,8,9},則[2,4] 位置上的值只有 3 種可能,為估算 S 的上界,取其最大的可能數,依次類推,於是有:
123 456 789 654 654 321 321 000 000
同樣地,考察第 2 宮和第 3 宮,會有:
123 456 789 654 654 321 321 321 321
繼續依次考察第 4 宮、第 2 列、第 3 列,會有:
123 456 789 654 654 321 321 321 321 665 000 000 543 000 000 421 000 000 333 000 000 222 000 000 111 000 000
繼續依次考察第 4 行、第 4 列、第 5 宮、第 5 行、第 6 行、第 5 列、第 6 列、第 7 列、第 7 行、第 8 行、第 9 行,會有:
123 456 789 654 654 321 321 321 321 665 654 321 543 543 321 421 421 321 333 333 321 222 222 221 111 111 111
這樣就得到了 S 的一個上界,即:
S < (6! · 6! · 6! · 3! · 3!) · (6! · 3! · 3! · 6!) · (5! · 4! · 3! · 3!) · (3! · 3! · 3! · 2! · 2!)
= 6!5 · 5! · 4! · 3!9 · 2!2 = 7205 · 120 · 24 · 69 · 4
= 193,491,763,200,000 · 116,095,057,920
= 22,463,437,455,746,924,544,000,000
於是
δ = 9! · S < 362,880 · 22,463,437,455,746,924,544,000,000
= 8,151,532,183,941,443,978,526,720,000,000
對比一下這些大數:
264 = 18,446,744,073,709,551,616
S < 22,463,437,455,746,924,544,000,000
296 = 79,228,162,514,264,337,593,543,950,336
δ < 8,151,532,183,941,443,978,526,720,000,000
可以看到 S 比 大很多數量級,因此上述的改進程式在m_soluSum 基礎上增加一個 m_bigSum 是不夠的,還需增加一個 m_bigbigSum(因為 S < 296);而如果直接對全空數獨題用程式去累計求 δ 的值,則還需要增加 m_bigbigbigSum(因為 δ > 296)。
具體估算一下程式對 1R 數獨題直接求 S 的值的時間。先看一下求一萬個解需要多長時間(對 1R 數獨題求解):
Run-run time: 2337 milliseconds; current solutions: 10000
steps: 0 # 0 # 445674
total solutions: 0 # 0 # 70003.
就用 2 秒一萬算,一小時 1800 萬,一天 4.32 億,一年是:
432,000,000 · 365 = 157,680,000,000
而 S ≈ 22,463,437,455,746,924,544,000,000
S 還是太大了。由上面的 ① 可以構造如下一個數獨題:
123 456 789 456 000 000 789 000 000 200 000 000 300 000 000 500 000 000 600 000 000 800 000 000 900 000 000
即第 1 行、第 1 宮、第 1列是滿填狀態,其餘都為空位,把這個數獨題稱作 1R1B1C 數獨題。假設這個題有 H 個解,由上述分析,有:
H = S / (6! · 6!) ≈ 43,332,248,178,524,160,000
而
264 = 18,446,744,073,709,551,616
296 = 79,228,162,514,264,337,593,543,950,336
雖然 H < 296,但是 2 秒一萬的求解速度還是太慢了,要求出 H 的具體值大約需要:
43,332,248,178,524,160,000 / 157,680,000,000 ≈ 274,811,315.19 年
約 2.75 億年啊!
這個問題還得尋求更好的數學方法來求解。
最後附上 SudokuSolver 2.5 的新增程式碼實現:
CQuizDealer 類聲明部分的修改
增加 runrun 介面:
void runrun(ulong newsum = 0);
新增步數和解數統計相關成員:
ulong m_soluSum; ulong m_bigSum; ulong m_bigbigSum; ulong m_steps; ulong m_bigSteps; ulong m_bigbigSteps;
incSteps 介面實現修改以及新增 incSolutions介面:
inline void incSteps() { if (m_steps != 0xFFFFFFFF) ++m_steps; else { m_steps = 0; if (m_bigSteps != 0xFFFFFFFF) printf("bigSteps:%u\n", ++m_bigSteps); else { m_bigSteps = 0; printf("bigbigSteps:%u\n", ++m_bigbigSteps); } } } inline void incSolutions() { ++s_soluSum; if (m_soluSum != 0xFFFFFFFF) ++m_soluSum; else { m_soluSum = 0; if (m_bigSum != 0xFFFFFFFF) printf("current big solution sum:%u\n", ++m_bigSum); else { m_bigSum = 0; printf("current big-big solution sum:%u\n", ++m_bigbigSum); } } }
其中 s_soluSum 是個靜態變數:
static ulong s_soluSum = 0;
相應介面實現小修改
parse 介面實現末尾的:
m_soluSum++;
換成
incSolutions();
adjust 介面實現內部也有一處同樣的替換。
runrun 介面實現
1 void CQuizDealer::runrun(ulong newsum) 2 { 3 if (m_state == STA_UNLOADED) { 4 printf("Quiz not loaded yet.\n"); 5 return; 6 } 7 s_soluSum = 0; 8 if (newsum == 0) 9 newsum = 0x10000; 10 clock_t begin = clock(); 11 while (m_state != STA_INVALID && s_soluSum < newsum) { 12 if (m_state == STA_LOADED) 13 parse(); 14 else if (m_state == STA_VALID) 15 adjust(); 16 else if (m_state == STA_DONE) { 17 if (m_stkSnap.empty()) 18 break; 19 m_state = STA_VALID; 20 nextGuess(); 21 } 22 } 23 std::cout << "Run-run time: " << clock() - begin << " milliseconds; current solutions: " << s_soluSum << std::endl; 24 std::cout << " steps: " << m_bigbigSteps << " # " << m_bigSteps << " # " << m_steps << std::endl; 25 std::cout << " total solutions: " << m_bigbigSum << " # " << m_bigSum << " # " << m_soluSum << ".\n"; 26 }
其他小修改
// 1.0 2021/9/20 // 2.0 2021/10/2 // 2.1 2021/10/4 // 2.2 2021/10/10 // 2.3 2021/10/17 // 2.4 2021/10/19 #define STR_VER "Sudoku Solver 2.5 2021/10/23 by readalps\n\n" void showOrderList() { ... printf("runrun <sum>: run till the end or certain new solutions met\n"); printf("bye: quit\n"); } void dealOrder(std::string& strOrder) { ...else if (matchPrefixEx(strOrder, "runrun ", strEx)) { ulong newsum = strtoul(strEx.c_str(), 0, 10); CQuizDealer::instance()->runrun(newsum); } else showOrderList(); }