朋友，這篇筆記看起來怎么樣呢？
前言：這是我 3月22日 寫的線性求解器 gecode 的筆記。我喜歡 gecode 的文檔，也很喜歡它簡潔流暢的 C++ 源碼。計劃趕不上變化快，當時興致勃勃的我，真的想不到，3月22日是我這五個月來最后一次寫這么工整的筆記;之后我便陷入一場極大的拉鋸戰中——無暇再打開這個筆記(甚至都給忘了，這幾天整理資料才想起)，但拉鋸戰最終收獲甚少;當時的我更想不到，現在提起「線性優化」這個課題，我的第一反應是拒絕。我很喜歡這個筆記，稍微閱讀我便可以把當時總結的重點全部回顧起來。這幾個月來的教訓是：不要認為自己可以什么都學、什么都涉獵，趕快定好一個細分方向，努力下去吧。接下來的很長一段時間，我可能都沒空更新這個筆記了。我當時的期望「gecode 是很優秀的線性求解器，但中文互聯網幾乎沒有其系統性教輔資料，我用研一下學期來寫一套」近期也不可能實現了。共勉。我會回來看你的，老伙計。
Gecode基礎知識
- Gecode基礎知識[1]
- SendMoreMoney[3]
- 必須注意的點[4]
- rel - 關系約束[5]
- distinct - 兩兩不同約束[6]
- branch - 定義分支規則[7]
- 關于 print[8]
- 以SendMoreMoney為例[2]
- 搜索引擎的使用[9]
- 異常處理[10]
- Gist - 可視化求解過程[11]
- 設定目標值[12]
- 筆記地址[13]
以SendMoreMoney為例
SendMoreMoney
我們給這八個字母賦值(賦一個數字)：S E N D M O R Y
使得等式成立：SEND + MORE = MONEY ，其中每個字母是不同的數字。
還不明白?告訴你答案，我們最后要把 SEND + MORE = MONEY 替換為 9567 + 1085 = 10652 。這就是求解成功了，其中 S 是 9 、 E 是 5 ...
必須注意的點
代碼在：https://gitee.com/piperliu/math_codes_economics_management/tree/master/gecode_MPG_notes/codes/03sendmoremoney/sendmoremoney.cpp[14]
- // 使用整數變量和整數約束
- // 我們需要 include gecode/int.hh 頭文件
- #include <gecode/int.hh>
- // 使用搜索引擎（求解邏輯），我們需要 gecode/search.hh 頭文件
- #include <gecode/search.hh>
- // 所有的 gecode 功能都在命名空間 Gecode 里
- // 所以，避免 Gecode::IntVarArray 這樣的麻煩
- // 我們 using namespace Gecode; 來簡化代碼
- using namespace Gecode;
- // 所有的模型都要繼承 Gecode::Space
- class SendMoreMoney : public Space {
- protected:
- // 這里聲明了一個整數變量數組 l
- IntVarArray l;
- public:
- // 在構造函數里聲明模型
- // l(*this, 8, 0, 9) 表示：
- // - 變量數組 l 是這個模型的（this）
- // - 變量數組 l 包含 8 個變量
- // - 這些變量的取值范圍是 [0, 9]
- SendMoreMoney(void) : l(*this, 8, 0, 9) {
- // 我們聲明了一些變量
- IntVar s(l[0]), e(l[1]), n(l[2]), d(l[3]),
- m(l[4]), o(l[5]), r(l[6]), y(l[7]);
- // 下面是一些約束的實現，先不用管具體含義
- // no leading zeros
- rel(*this, s, IRT_NQ, 0);
- rel(*this, m, IRT_NQ, 0);
- // all letters distinct
- distinct(*this, l);
- // linear equation
- // 這里是一些 Args 變量
- // 我把 Var 理解為決策變量，或者自變量
- // 而 Args 是根據 Var 變動的，中介變量
- // 本身并不主動變化，受 Var 們變動影響
- IntArgs c(4+4+5); IntVarArgs x(4+4+5);
- c[0]=1000; c[1]=100; c[2]=10; c[3]=1;
- x[0]=s; x[1]=e; x[2]=n; x[3]=d;
- c[4]=1000; c[5]=100; c[6]=10; c[7]=1;
- x[4]=m; x[5]=o; x[6]=r; x[7]=e;
- c[8]=-10000; c[9]=-1000; c[10]=-100; c[11]=-10; c[12]=-1;
- x[8]=m; x[9]=o; x[10]=n; x[11]=e; x[12]=y;
- linear(*this, c, x, IRT_EQ, 0);
- // 定義了分支定界的規則
- // post branching
- branch(*this, l, INT_VAR_SIZE_MIN(), INT_VAL_MIN());
- }
- // 必須聲明一個 copy constructor
- // 其中一定一定要有對于變量的更新
- // search support
- SendMoreMoney(SendMoreMoney& s) : Space(s) {
- l.update(*this, s.l);
- }
- // 以及聲明一個 copy 方法（用于search）
- virtual Space* copy(void) {
- return new SendMoreMoney(*this);
- }
- // print solution
- void print(void) const {
- std::cout << l << std::endl;
- }
- };
- // main function
- int main(int argc, char* argv[]) {
- // create model and search engine
- SendMoreMoney* m = new SendMoreMoney;
- DFS<SendMoreMoney> e(m);
- delete m;
- // search and print all solutions
- while (SendMoreMoney* s = e.next()) {
- s->print(); delete s;
- }
- return 0;
- }
如上，我們實現了一個問題，我在中文注釋中做出了解釋。值得注意的有：
- 所有的 gecode 功能都在命名空間 Gecode 里
- 所有的模型都要繼承 Gecode::Space
- 除了基本的構造函數，必須聲明一個 copy constructor 、以及聲明一個 copy 方法(用于search)
- 我們必須在 copy constructor 對決策變量進行更新(比如 l.update(*this, s.l);)，這尤為重要，否則程序不會報錯，而求解過程是錯誤的。原因在于，我們的求解器在求解過程中，會不斷賦值 Space 實例，如果不在 copy constructor 中更新決策變量，則無法進行回溯、保存等操作
我把 Var 理解為決策變量，或者自變量;而 Args 是根據 Var 變動的，中介變量，本身并不主動變化，受 Var 們變動影響。
注意到無論是變量還是約束，其聲明時，第一個參數都 *this 。你可以如此理解：變量如小孩子，小孩子出生，必須明確 ta 的監護人是誰，否則，接下來小孩子沒有依靠，這個小孩子也就活不下去。在這里，我們的變量和約束都依靠其存在的問題中，即類 SendMoreMoney 。
rel - 關系約束
- rel(*this, s, IRT_NQ, 0);
上述語句表示：變量 s 不等于 0，其中：
- IRT_NQ 表示不等于，我們將在后文見到所有的邏輯符號表達
- rel 是 relation 的意思
distinct - 兩兩不同約束
- distinct(*this, l);
我們這里不具體闡述 distinct 的用法，只說明其作用：
- l 數組中變量要兩兩不同(pairwise different)
其他約束不再詳述，以后會提到。
branch - 定義分支規則
在求解過程中，不同的搜索方式對于求解效果有著天差地別的影響。因此定義合適的 branch 分支規則尤為重要。
- branch(*this, l, INT_VAR_SIZE_MIN(), INT_VAL_MIN());
這里表示，對于變量數組 l 中的變量，我們首先選擇對其中搜索范圍最小的變量進行搜索(INT_VAR_SIZE_MIN())，此外如果選擇了一個變量，我們將現有的最小值賦給它(INT_VAL_MIN())。
其他規則以及關于 branch 的用法將在后文討論。
關于 print
在求解后， gecode 會自動調用 print 輸出。我們可以自由發揮，重載 print 函數，讓解是我想要的形式。
- // print solution
- void print(void) const {
- std::cout << l << std::endl;
- }
求解后的效果：
- {9, 5, 6, 7, 1, 0, 8, 2}
表示 S 是 9 、 E 是 5 ...
搜索引擎的使用
我們在上面定義了我們的問題(包含變量與約束)，我們在程序中稱之為一個“模型model”。
這個模型對象中存儲了當前變量的范圍，我將其理解為搜索中的一些節點的集合(在一些運算后產生的更小的、更精準的解空間)。
那么如何實現一顆搜索樹呢?在 main 方法中，我們將其實現。
- // main function
- int main(int argc, char* argv[]) {
- // create model and search engine
- SendMoreMoney* m = new SendMoreMoney;
- // 第一個 print()
- std::cout << "1" << std::endl;
- m->print();
- m->status();
- // 第二個 print()
- std::cout << "2" << std::endl;
- m->print();
- DFS<SendMoreMoney> e(m);
- // 第三個 print()
- std::cout << "3" << std::endl;
- m->print();
- delete m;
- // search and print all solutions
- while (SendMoreMoney* s = e.next()) {
- // 第四個 print()
- std::cout << "4" << std::endl;
- s->print(); delete s;
- }
- return 0;
- }
你可以看到，我在上面安插了 4 個 print，來看看當前我們搜索到的范圍。輸入如下。
- 1
- {[1..9], [0..9], [0..9], [0..9], [1..9], [0..9], [0..9], [0..9]}
- 2
- {9, [4..7], [5..8], [2..8], 1, 0, [2..8], [2..8]}
- 3
- {9, [4..7], [5..8], [2..8], 1, 0, [2..8], [2..8]}
- 4
- {9, 5, 6, 7, 1, 0, 8, 2}
當我們聲明這個問題時(初始化模型)，我們的模型根據約束做了一些最基本的運算，在這里體現為：
- S的范圍應該是[1, 9]
- M的范圍應該是[1, 9]
這很好理解，因為我們有約束 rel(*this, s, IRT_NQ, 0);， rel(*this, m, IRT_NQ, 0); 。
只有當我們調用 status() 后，約束傳播constraint propagation才會起作用。
status函數的注釋原文為 Query space status Propagates the space until fixpoint or failure; updates the statistics information stat; and:
- if the space is failed, SpaceStatus::SS_FAILED is returned.
- if the space is not failed but the space has no brancher left, SpaceStatus::SS_SOLVED is returned.
- otherwise, SpaceStatus::SS_BRANCH is returned.
我的理解是，對現有的解空間做一次最基本的探索，直到需要 branch 做決策的點或者探索失敗。
值得注意的是，如果我們取消上述 m->status(); 這句話，那么 3 對應的 print() 依然是 {9, [4..7], [5..8], [2..8], 1, 0, [2..8], [2..8]} 。這說明為我們的模型聲明搜索引擎時，會調用一次 status() 。
如下，在 while 中，我們希望返回每一次探索的解。值得注意的是，解也是一個模型(問題實例)。
- while (SendMoreMoney* s = e.next()) {
- // 第四個 print()
- std::cout << "4" << std::endl;
- s->print(); delete s;
- }
在這里，e是基于DFS深度優先搜索的搜索引擎，e.next()則是基于DFS規則做出的下一個探索方向的探索結果。如果探索到頭了，那么e.next()會返回NULL，此時 while 停止，跳出循環。
異常處理
實際上，為了捕捉 Gecode 的異常，使用 try {} catch () {} 是一個很好的選擇。
- int main(int argc, char* argv[]) {
- try {
- // main 中的內容
- } catch (Exception e) {
- std::cerr << "Gecode exception: " << e.what() << std::endl;
- return 1;
- }
- return 0;
- }
注意，這里的 Exception 是 Gecode 命名空間中的異常Gecode::Exception。
不管是本筆記還是教程文件MPG，都沒有使用 try catch (為了易讀性)，但是在時間使用中，還是應該具有如此的良好習慣。
Gist - 可視化求解過程
程序見https://gitee.com/piperliu/math_codes_economics_management/tree/master/gecode_MPG_notes/codes/03sendmoremoney/sendmoremoneyGist_1.cpp[15]
- #include <gecode/gist.hh>
- using namespace Gecode;
- ...
- int main(int argc, char* argv[]) {
- SendMoreMoney* m = new SendMoreMoney;
- Gist::dfs(m);
- delete m;
- return 0;
- }
我們引入 gist.hh 頭文件，在 main 中較為簡單地可視化 dfs 過程。
可以用鼠標探索
關于 gist 功能，還有其他調用方法。gist功能較為強大，這里不再贅述。
此外，還可用手動增加現實解的功能，代碼見https://gitee.com/piperliu/math_codes_economics_management/tree/master/gecode_MPG_notes/codes/03sendmoremoney/sendmoremoneyGist_2.cpp[16]
添加了 inspect
設定目標值
設定目標值時，我們需要聲明 cost 函數以及用變量表示搜索過程中變量的傳遞(用于目標函數的計算)。這里不詳細討論。
這并不復雜，我們將在之后的筆記中專門討論。
