戰(zhàn)：從集合語法到求解器選型與靈敏度分析)
簡(jiǎn)介L(zhǎng)INGO線性與非線性全局優(yōu)化器是一款功能強(qiáng)大的數(shù)學(xué)建模優(yōu)化軟件本資源提供其破解安裝版本面向工程、經(jīng)濟(jì)、金融、運(yùn)籌學(xué)等領(lǐng)域的研究者、工程師及建模競(jìng)賽參與者用于求解線性、非線性與整數(shù)規(guī)劃問題。軟件內(nèi)置專門的建模語言用戶可通過簡(jiǎn)潔語句定義變量、目標(biāo)函數(shù)和約束快速搭建大規(guī)模優(yōu)化模型求解器性能強(qiáng)勁支持混合整數(shù)規(guī)劃和結(jié)果敏感性分析還提供數(shù)據(jù)接口與報(bào)告生成功能便于多場(chǎng)景建模與對(duì)比。壓縮包為rar格式大小約15.9MB輕量易獲取目前已有650人學(xué)習(xí)/下載適合從入門到進(jìn)階的不同層級(jí)用戶。借助此資源可直接體驗(yàn)LINGO 11版本的增強(qiáng)求解能力與界面優(yōu)化結(jié)合官方文檔和實(shí)例在供應(yīng)鏈調(diào)度、投資組合優(yōu)化、能源系統(tǒng)設(shè)計(jì)等實(shí)際場(chǎng)景中快速構(gòu)建最優(yōu)解方案提升建模與決策效率。1. LINGO 不是計(jì)算器它到底在解什么LINGO 這個(gè)工具表面看是把數(shù)學(xué)公式敲進(jìn)去、點(diǎn) Solve 就出結(jié)果的“計(jì)算器”但真正拉開建模差距的是它對(duì)集合、數(shù)據(jù)段和求解器分組這三個(gè)設(shè)計(jì)點(diǎn)的處理。同一道排產(chǎn)題有人用幾十行腳本勉強(qiáng)跑通有人用幾行集合寫法在秒級(jí)拿到最優(yōu)解差別不在數(shù)學(xué)功底而在建模方式。這篇文章不打算按教科書順序講語法而是從建模語法、求解器選型、報(bào)告解讀到批量聯(lián)調(diào)完整走一遍我在真實(shí)項(xiàng)目里驗(yàn)證過的流程。適合手頭有優(yōu)化問題的從業(yè)者——供應(yīng)鏈調(diào)度、生產(chǎn)計(jì)劃、定價(jià)決策都算——也適合想快速驗(yàn)證算法思路的同學(xué)。建議開著 LINGO 邊看邊跑大部分示例改兩行數(shù)據(jù)就能直接復(fù)現(xiàn)。2. 集合與數(shù)據(jù)段把模型從紙面搬進(jìn) LINGO 的兩種寫法2.1 原始集與派生集先把“對(duì)誰建?！睂懬宄﨤INGO 和直接寫 LP/MPS 文件的習(xí)慣很不一樣。從紙面公式搬到 LINGO第一件事不是寫目標(biāo)函數(shù)而是定義集合。集合相當(dāng)于告訴求解器我要對(duì)哪些對(duì)象建模型每個(gè)對(duì)象帶哪些屬性。以最經(jīng)典的運(yùn)輸問題為例需求端是 4 個(gè)客戶供給端是 3 個(gè)倉庫每個(gè)倉庫有供應(yīng)量每個(gè)客戶有需求量每條倉庫到客戶的路徑有單位運(yùn)輸成本和運(yùn)輸量。MODEL: SETS: WAREHOUSE /W1, W2, W3/: SUPPLY; CUSTOMER /C1, C2, C3, C4/: DEMAND; ROUTE(WAREHOUSE, CUSTOMER): COST, VOLUME; ENDSETS DATA: SUPPLY 30 40 50; DEMAND 20 25 30 45; COST 2 4 5 3 3 5 4 6 4 2 3 5; ENDDATA MIN SUM(ROUTE(I, J): COST(I, J) * VOLUME(I, J)); FOR(CUSTOMER(J): SUM(WAREHOUSE(I): VOLUME(I, J)) DEMAND(J)); FOR(WAREHOUSE(I): SUM(CUSTOMER(J): VOLUME(I, J)) SUPPLY(I)); FOR(ROUTE(I, J): VOLUME(I, J) 0); END這段代碼里的三類寫法建議記牢。第一原始集定義用名字 /成員/ : 屬性成員之間用逗號(hào)或空格隔開屬性跟在冒號(hào)后面一個(gè)屬性就是一個(gè)和集合同長(zhǎng)度的參數(shù)或變量數(shù)組。第二派生集ROUTE(WAREHOUSE, CUSTOMER)表示兩個(gè)原始集的笛卡爾積展開后共 12 個(gè)成員COST和VOLUME都掛在這個(gè)派生集上分別讀入 12 個(gè)成本數(shù)據(jù)和生成 12 個(gè)決策變量。第三FOR和SUM是 LINGO 的循環(huán)與求和函數(shù)FOR(CUSTOMER(J): ...)表示對(duì)每個(gè)客戶單獨(dú)生成一條等式約束SUM(WAREHOUSE(I): ...)表示對(duì)倉庫下標(biāo)求和。我一般會(huì)先在草稿紙上畫一張“對(duì)象關(guān)系圖”把集合畫成方框、屬性畫成方框里的行、跨集合的連線就是派生集。這樣做直接的好處是寫FOR時(shí)不容易搞混下標(biāo)。比如上面第二個(gè)FOR里外層循環(huán)下標(biāo)是I倉庫內(nèi)部求和下標(biāo)是J客戶意思是每個(gè)倉庫發(fā)出去的總量等于它的供應(yīng)量如果把兩者位置對(duì)調(diào)模型含義就從“運(yùn)出量守恒”變成“運(yùn)入量等于別的東西”LINGO 不會(huì)報(bào)警告只有結(jié)果跟業(yè)務(wù)對(duì)不上時(shí)你才會(huì)發(fā)現(xiàn)。2.2 數(shù)據(jù)段與外部文件把參數(shù)和模型結(jié)構(gòu)分開上一節(jié)的運(yùn)輸問題里數(shù)據(jù)全寫在DATA/ENDDATA段中這是我最推薦的起始方式模型結(jié)構(gòu)一眼能看完改數(shù)據(jù)也方便。但數(shù)據(jù)量超過幾十個(gè)、或者同一套模型要套用到多個(gè)區(qū)域時(shí)就應(yīng)該把參數(shù)挪到外部文件。LINGO 的FILE函數(shù)就是干這個(gè)的。DATA: SUPPLY FILE(sup_data.txt); DEMAND FILE(dem_data.txt); COST FILE(cost_data.txt); ENDDATAsup_data.txt里只需要按集合成員順序?qū)憯?shù)值例如30 40 50文件里可以有換行和空格數(shù)值之間不需要逗號(hào)。這樣模型文件里只剩集合和約束結(jié)構(gòu)數(shù)據(jù)歸數(shù)據(jù)。我做過的一次模擬項(xiàng)目 X就是同一套庫存補(bǔ)貨模型把華東、華南、華北三組數(shù)據(jù)分別存成三個(gè)文本跑三次模型一行沒改。把數(shù)據(jù)外置還有一個(gè)隱藏優(yōu)點(diǎn)做參數(shù)掃描時(shí)你只需要用腳本替換文本里的某個(gè)值然后重新求解。后面第 6 章我會(huì)給完整的批量做法。如果你想讓結(jié)果導(dǎo)出成可讀文件DATA段里也可以寫TEXT(result.txt)把某個(gè)集合屬性的值按成員順序?qū)懙轿谋疚募?。這個(gè)函數(shù)我在前面的選址問題里常用至少比從密密麻麻的 Solution Report 里手工復(fù)制要穩(wěn)。2.3 從紙面公式到 LINGO 表達(dá)的三個(gè)習(xí)慣最后補(bǔ)三個(gè)我踩過坑練出來的習(xí)慣。第一目標(biāo)函數(shù)寫在最前面MIN或MAX開頭接下來所有等號(hào)、不等式、BND都是約束。第二約束盡量加行標(biāo)簽比如FOR(CUSTOMER(J): [DEMAND_BALANCE] SUM(...) DEMAND(J));加標(biāo)簽之后報(bào)告里的 Row 編號(hào)就會(huì)帶著名字排查結(jié)果時(shí)能少一半時(shí)間去對(duì)行號(hào)。第三決策變量和參數(shù)不要混在一個(gè)集合屬性里參數(shù)在DATA段賦值變量只聲明不賦值如果非要給變量初值可以在DATA段賦但別把初值當(dāng)業(yè)務(wù)參數(shù)用否則下次換數(shù)據(jù)時(shí)容易漏改。3. 求解器選型LP、MIP、NLP 跑不對(duì)時(shí)先查這三處3.1 線性模型與整數(shù)模型變量類型聲明和 Gap 設(shè)置LINGO 會(huì)根據(jù)模型自動(dòng)分派求解器純線性連續(xù)問題是 LP帶整數(shù)變量是 MIP目標(biāo)或約束里有非線性表達(dá)式就是 NLP。自動(dòng)分派省事但很多“跑出來不對(duì)”的翻車現(xiàn)場(chǎng)恰恰是因?yàn)槟Ｐ捅蛔詣?dòng)分派到了錯(cuò)誤的求解器。最典型的是整數(shù)變量沒聲明模型被當(dāng)成 LP 解得到小數(shù)方案。設(shè)施選址問題能同時(shí)展示 LP 到 MIP 的差別。MODEL: SETS: PLANT /P1, P2, P3/: FIXED, CAPACITY, OPEN; CUSTOMER /C1, C2, C3, C4/: DEMAND; ROUTE(PLANT, CUSTOMER): TRANSCOST, SHIP; ENDSETS DATA: FIXED 1000 800 1200; CAPACITY 500 400 600; DEMAND 150 180 120 200; TRANSCOST 5 4 6 3 7 3 5 6 4 5 4 7; ENDDATA MIN SUM(PLANT(I): FIXED(I) * OPEN(I)) SUM(ROUTE(I, J): TRANSCOST(I, J) * SHIP(I, J)); FOR(CUSTOMER(J): SUM(PLANT(I): SHIP(I, J)) DEMAND(J)); FOR(PLANT(I): SUM(CUSTOMER(J): SHIP(I, J)) CAPACITY(I) * OPEN(I)); FOR(ROUTE(I, J): SHIP(I, J) 0); FOR(PLANT(I): BIN(OPEN(I))); END關(guān)鍵在最后兩行。BIN(OPEN(I))把OPEN聲明成 0/1 變量模型才是 MIP倒數(shù)第二條容量約束 CAPACITY(I) * OPEN(I)的含義是OPEN為 1 時(shí)該廠產(chǎn)能正常生效為 0 時(shí)上限變成 0運(yùn)輸量被迫為 0。這個(gè)“乘法開關(guān)”模式幾乎適用于所有帶固定成本的選擇問題。MIP 跑得慢或者太久不出最優(yōu)解時(shí)先別急著加機(jī)器去 Options 里的 MIP Solver 參數(shù)看 Dual Gap。Gap 表示當(dāng)前最好的整數(shù)可行解與松弛上界之間的相對(duì)差距設(shè)成 0 表示一定要證明最優(yōu)才停代價(jià)是分支定界可能跑幾個(gè)小時(shí)設(shè)成 1% 通常能在十幾分鐘內(nèi)給出工程上可接受的解。我在做某跨平臺(tái)系統(tǒng)的排產(chǎn)模塊時(shí)就把 Gap 設(shè)成 0.5%業(yè)務(wù)方完全接受求解時(shí)間從 40 分鐘壓到 3 分鐘。要注意如果模型本身只有幾千個(gè)變量Gap 設(shè)大點(diǎn)無所謂一旦變量上十萬Gap 開 1% 和開 5% 的耗時(shí)差能到 10 倍以上量力設(shè)。3.2 非線性模型與全局求解器默認(rèn)只是局部收斂非線性模型的默認(rèn)行為容易讓人誤判。LINGO 碰到非線性目標(biāo)或約束時(shí)默認(rèn)走 NLP 求解器本質(zhì)是局部搜索從一個(gè)初始點(diǎn)出發(fā)沿可行方向逐步改進(jìn)找到一個(gè)“站在附近看確實(shí)最優(yōu)”的解就停。如果模型是凸的局部即全局但現(xiàn)實(shí)中定價(jià)、選址、排程模型很少有凸的于是局部解就變成一個(gè)大坑。一個(gè)能快速?gòu)?fù)現(xiàn)的示例是圓柱儲(chǔ)罐用料最小化給定容積 100 立方米求半徑 R 和高度 H 使表面積最小。MODEL: ! 目標(biāo)最小化圓柱儲(chǔ)罐的用料表面積 V 100; ! 容量固定 100 立方米 MIN 2 * 3.14159 * R * R 2 * 3.14159 * R * H; BND(1, R, 10); BND(1, H, 50); 3.14159 * R * R * H V; END這個(gè)模型手算就能驗(yàn)證把 H 換成 100 / (3.14159 * R^2)目標(biāo)函數(shù)變成一元函數(shù)最優(yōu)在 R 約 2.5、H 約 5.0材料面積大概 118 平方米。默認(rèn) NLP 從初值 R1、H20 也能找到這個(gè)點(diǎn)因?yàn)榭尚杏蚯嫔线@個(gè)點(diǎn)是唯一谷底。但如果把目標(biāo)函數(shù)換成帶SIN或EXP的生產(chǎn)計(jì)劃模型谷底可能就不止一個(gè)。這時(shí)候要在 LINGO 菜單里打開 Options→Global Solver勾選啟用全局求解器它會(huì)在后臺(tái)做類似分支定界的過程反復(fù)逼近并證明全局最優(yōu)。代價(jià)是時(shí)間變量一多、函數(shù)一復(fù)雜跑幾十分鐘甚至幾天都正常。我一般先把局部解跑出來看結(jié)果是否在業(yè)務(wù)合理區(qū)間再用全局求解器驗(yàn)證一次而不是一開始就全局硬算。3.3 三個(gè)參數(shù)與數(shù)值尺度TIMLIM、Scale Model 和線性化選項(xiàng)求解器參數(shù)里我最常動(dòng)的是三個(gè)。第一個(gè)是時(shí)間上限 TIMLIM單位秒設(shè)到 300 到 600 很常見它防止 MIP 或全局求解器無限跑下去。第二個(gè)是迭代上限 ITERLIM只對(duì)單純形法有意義內(nèi)點(diǎn)法和分支定界主要看時(shí)間。第三個(gè)是 Scale Model 選項(xiàng)遇到目標(biāo)上百萬、約束系數(shù)只有 0.001 這類數(shù)值尺度失衡時(shí)勾上它經(jīng)常會(huì)讓 Infeasibilities 從 0.01 突然變成 0.000000。原理很簡(jiǎn)單矩陣?yán)锿瑫r(shí)存在 1e6 和 1e-4 級(jí)的數(shù)求解器的內(nèi)積、求逆精度會(huì)吃緊縮放后條件數(shù)改善結(jié)果立刻穩(wěn)。另外 LINGO 有自動(dòng)線性化能力小規(guī)模非線性表達(dá)式會(huì)在預(yù)求解階段嘗試線性化但一旦變量超過幾百個(gè)就別指望它全自動(dòng)人工引入輔助變量改寫反而更快這部分在第 5 章會(huì)給改寫例子。4. 報(bào)告與靈敏度分析從 Solution Report 里讀出真實(shí)邊界4.1 先看懂求解狀態(tài)和三個(gè)關(guān)鍵字段每次求解完彈出的 Solution Report大多數(shù)人只掃一眼目標(biāo)值就關(guān)掉。我建議至少看四個(gè)地方。首先是報(bào)告第一行Global optimal solution found說明已經(jīng)證明全局最優(yōu)Local optimum found說明只是局部解Feasible solution found更弱只保證找到可行解。這三者差一個(gè)級(jí)別業(yè)務(wù)決策的置信度就差一個(gè)檔次。其次看Infeasibilities它不為 0 說明模型找不到嚴(yán)格滿足所有約束的點(diǎn)這時(shí)目標(biāo)值再好看都不能用。下面是一個(gè)簡(jiǎn)化后的 LP 小案例報(bào)告用來解釋字段含義不必對(duì)應(yīng)前面的某個(gè)具體模型。Global optimal solution found. Objective value: 870.0000 Objective bound: 870.0000 Infeasibilities: 0.000000 Total solver iterations: 5 Variable Value Reduced Cost X1 40.00000 0.000000 X2 20.00000 5.000000 Row Slack or Surplus Dual Price 1 870.0000 1.000000 2 0.000000 7.000000 3 10.00000 0.000000變量區(qū)里的Reduced Cost對(duì)基變量通常為 0對(duì)非基變量表示如果把這個(gè)變量強(qiáng)行加進(jìn)最優(yōu)解一個(gè)單位目標(biāo)值會(huì)惡化多少。Row 區(qū)里的Slack or Surplus是約束的剩余量等于 0 說明該約束緊貼邊界是真正的限制因素大于 0 說明約束還有余量暫時(shí)不影響目標(biāo)。Dual Price就是影子價(jià)格下面單獨(dú)說。Row 后面的編號(hào)對(duì)應(yīng)模型里的約束順序所以我在第 2.3 節(jié)強(qiáng)調(diào)要加行標(biāo)簽就是為了讓這里的編號(hào)能一眼對(duì)上約束名。我在實(shí)戰(zhàn)里的習(xí)慣是改動(dòng)模型后先看求解狀態(tài)那行字再看 Infeasibilities最后才看目標(biāo)值。狀態(tài)不對(duì)時(shí)數(shù)值再漂亮也只是黑匣子輸出。4.2 靈敏度 Range 表影子價(jià)格能走多遠(yuǎn)默認(rèn)報(bào)告不會(huì)給完整的靈敏度 Range 表。要在求解前把 Options 里 General Solver 的 Dual Computations 設(shè)成 Prices Ranges重新求解后才會(huì)出現(xiàn)兩組區(qū)間目標(biāo)系數(shù)變化范圍和右側(cè) RHS 變化范圍。目標(biāo)系數(shù)范圍表長(zhǎng)這樣。變量當(dāng)前系數(shù)允許增加允許減少X150.0010.0020.00X240.008.006.00RHS 范圍表單獨(dú)一張。約束當(dāng)前 RHS允許增加允許減少約束2100.005.008.00約束380.0012.0020.00這兩張表回答的問題是我不重新求解的情況下業(yè)務(wù)參數(shù)在多大范圍內(nèi)波動(dòng)當(dāng)前這組最優(yōu)決策仍然是“結(jié)構(gòu)最優(yōu)”。注意“結(jié)構(gòu)最優(yōu)”不等于“變量值不變”它只表示最優(yōu)基的變量集合不變。一旦某個(gè)系數(shù)的實(shí)際增幅超出允許增加最優(yōu)基會(huì)變影子價(jià)格會(huì)跳到另一個(gè)值這時(shí)候再去引用原 Dual Price 就作廢了。影子價(jià)格業(yè)務(wù)化是最容易翻車的地方。某一次我給一個(gè)倉儲(chǔ)項(xiàng)目做資源評(píng)估報(bào)告里約束 2 的 Dual Price 是 7意思是該約束右側(cè)再給 1 個(gè)單位目標(biāo)值能改善 7。業(yè)務(wù)方很興奮說多放 10 個(gè)貨架就多賺 70 萬。我拿 Range 表一看這個(gè)約束的 RHS 允許增加只有 5第 6 個(gè)貨架開始 Shadow Price 已經(jīng)不是 7 了。所以我把影子價(jià)格當(dāng)“邊際價(jià)值信號(hào)”不當(dāng)“線性外推的乘法口訣表”并且只對(duì)它生效的范圍負(fù)責(zé)。4.3 用手工可解案例做回歸報(bào)告說最優(yōu)不等于邏輯對(duì)還有一層驗(yàn)證是很多人跳過的報(bào)告說最優(yōu)不等于模型建得沒問題。我習(xí)慣在正式跑大規(guī)模數(shù)據(jù)前先構(gòu)造一個(gè)規(guī)模小到能手工算的案例塞進(jìn)同一套模型確認(rèn)結(jié)果和手算一致。比如運(yùn)輸問題我做過 2 個(gè)倉庫送 2 個(gè)客戶的小案例手算最優(yōu)值再讓模型跑對(duì)上了才換真實(shí)數(shù)據(jù)。這件事花 10 分鐘能擋住一大半“集合寫錯(cuò)但報(bào)告照出最優(yōu)解”的隱性 bug。5. 避坑實(shí)錄五個(gè)反復(fù)出現(xiàn)的建模翻車場(chǎng)景與處置5.1 整數(shù)聲明缺失解出 13.7 臺(tái)設(shè)備業(yè)務(wù)當(dāng)場(chǎng)不認(rèn)現(xiàn)象生產(chǎn)計(jì)劃模型跑完報(bào)告里生產(chǎn)量是 13.7 臺(tái)、開啟 0.4 條產(chǎn)線業(yè)務(wù)方直接拒收結(jié)果。原因LINGO 變量默認(rèn)是連續(xù)變量GIN、BIN沒寫模型被當(dāng)成 LP 解了。解決凡是決策量必須是整數(shù)的聲明時(shí)給集合屬性補(bǔ)一行GIN(PRODUCE(P))或BIN(OPEN(I))。跑完先看報(bào)告底部的 Model Statistics核對(duì)整數(shù)變量個(gè)數(shù)。我后來養(yǎng)成的習(xí)慣是需求評(píng)審階段就把“哪些量必須是整數(shù)”寫進(jìn)建模文檔避免寫代碼時(shí)漏掉。5.2 集合維度錯(cuò)配Error 72 的鍋十有八九在這現(xiàn)象加載模型報(bào) Runtime Error錯(cuò)誤碼 72指向某一行FOR或SUM。原因下標(biāo)集合寫錯(cuò)了維度。比如在FOR(CUSTOMER(J): ...)內(nèi)部卻用了SUM(WAREHOUSE(J): SHIP(I, J))把客戶集合的成員J當(dāng)成了倉庫集合的成員用。解決先定位報(bào)錯(cuò)行檢查每個(gè)循環(huán)變量所屬的集合。FOR(A(I): ...)里可以安全使用的下標(biāo)只有I和它所在約束里通過SUM跨集合引用的其他集合下標(biāo)。一個(gè)笨但有效的辦法把模型里所有FOR和SUM的下標(biāo)變量名列出來逐個(gè)對(duì)照 SETS 定義。另外派生集下標(biāo)的順序要和定義時(shí)一致ROUTE(WAREHOUSE, CUSTOMER)里第一維是倉庫第二維是客戶寫ROUTE(J, I)卻把 I 當(dāng)倉庫、J 當(dāng)客戶就會(huì)維度錯(cuò)配。5.3 IF 和 SMAX 把問題拆成非光滑懸崖現(xiàn)象模型能解但換一組合理初值結(jié)果跳到一個(gè)完全不同的解業(yè)務(wù)怎么調(diào)參數(shù)都達(dá)不到報(bào)告里的“最優(yōu)”。原因目標(biāo)或約束里用了IF、SMAX、ABS這類函數(shù)可行域被切成不連續(xù)的多塊非凸不可導(dǎo)默認(rèn) NLP 局部搜索只能爬到其中一個(gè)谷底。解決把非光滑項(xiàng)改寫成線性不等式。比如約束里出現(xiàn)SMAX(0, PRICE(I) - BID(I))可以引入輔助變量OVER(I)替換成下面兩條約束。OVER(I) PRICE(I) - BID(I); OVER(I) 0;原來用SMAX的地方直接換成OVER(I)。數(shù)學(xué)上等價(jià)數(shù)值上卻從“不可導(dǎo)懸崖”變成“兩條線性邊界”全局求解器才真正能用。類似地IF(COND, A, B)可以分解成一個(gè) 0/1 變量加兩個(gè)大 M 不等式。需要提醒的是大 M 不要拍腦袋給 1e6取“業(yè)務(wù)上該量級(jí)最大可能值的 3 倍”就夠太大反而把數(shù)值條件數(shù)拖壞。5.4 數(shù)值尺度懸殊手算滿足約束報(bào)告卻報(bào) Infeasibilities現(xiàn)象目標(biāo)函數(shù)動(dòng)輒上百萬約束系數(shù)小到 0.0001報(bào)告里 Infeasibilities 顯示 0.01可人工逐個(gè)約束驗(yàn)算發(fā)現(xiàn)都滿足。原因模型矩陣的條件數(shù)太差求解器在數(shù)值容差內(nèi)判斷違反約束產(chǎn)生了“數(shù)值偽不可行”。這在小數(shù)單位和萬元單位混用的成本核算模型里很常見。解決先統(tǒng)一業(yè)務(wù)單位把元換成萬元或把千克換成噸讓目標(biāo)值和系數(shù)都落在 1e-2 到 1e4 之間再去 Options 里勾選 Scale Model。我在某物流成本模型里遇到過 Infeasibilities0.03把單箱成本從“元”改成“萬元”后直接歸零一秒沒多算。如果改完還不行重點(diǎn)查是否有同時(shí)除以兩個(gè)量級(jí)相差很大的分母的表達(dá)式。5.5 局部最優(yōu)當(dāng)全局最優(yōu)同一個(gè)模型換初值換結(jié)果現(xiàn)象同一個(gè) NLP 模型同事用另一套求解器算出來的目標(biāo)值明顯比我的好換初始點(diǎn)再跑結(jié)果又變一個(gè)樣。原因默認(rèn) NLP 是局部搜索非凸問題從不同起點(diǎn)收斂到不同局部解根本沒有“唯一答案”。解決重要決策用全局求解器跑一遍同時(shí)設(shè)時(shí)間上限比如 600 秒看它找到的目標(biāo)值和默認(rèn) NLP 是否一致。如果全局求解時(shí)間完全不可接受退一步做多起點(diǎn)掃描人工挑 5 到 10 個(gè)有業(yè)務(wù)含義的初始點(diǎn)跑局部 NLP記錄每個(gè)起點(diǎn)對(duì)應(yīng)的目標(biāo)值和決策量取業(yè)務(wù)上最合理的那組并在文檔里注明這是“驗(yàn)證過的局部解”不是全局最優(yōu)。從那以后我對(duì)“最優(yōu)”這個(gè)詞都加限定語是全局證明還是局部搜索加經(jīng)驗(yàn)判斷。5.6 排查速查表癥狀最先查的語法/參數(shù)本章對(duì)應(yīng)場(chǎng)景解出小數(shù)GIN / BIN 聲明5.1Runtime Error 72FOR/SUM 的下標(biāo)集合匹配5.2換初值結(jié)果大變IF/SMAX/ABS 非線性改寫5.3Infeasibilities 不歸零Scale Model、單位歸一化5.4報(bào)告只有 Local optimum開 Global Solver5.56. 把 LINGO 嵌進(jìn)工作流Excel 聯(lián)調(diào)與批量參數(shù)掃描6.1 用 OLE 讓 Excel 當(dāng)數(shù)據(jù)界面Excel 里先給區(qū)域命名LINGO 的OLE就能按名字讀寫。數(shù)據(jù)段讀入?yún)?shù)CALC 段寫回結(jié)果業(yè)務(wù)方改 Excel 就行不碰模型文件。DATA: SUPPLY OLE(C:\supply_plan.xlsx, supply); DEMAND OLE(C:\supply_plan.xlsx, demand); COST OLE(C:\supply_plan.xlsx, cost); ENDDATA CALC: OLE(C:\supply_plan.xlsx, ship_out) SHIP; ENDCALC注意路徑里別用中文名Excel 文件要事先打開或在名稱管理器里確認(rèn)區(qū)域存在。第一次聯(lián)調(diào)失敗九成是 Excel 里的命名區(qū)域尺寸和集合成員數(shù)不一致。提示Excel 命名區(qū)域的尺寸必須和集合成員數(shù)一致否則 OLE 讀入會(huì)在運(yùn)行時(shí)報(bào)告維度不匹配先把區(qū)域重新框選再命名一次。6.2 參數(shù)掃描改文件、跑求解、讀結(jié)果需要做敏感性掃描時(shí)我不喜歡手動(dòng)改 LINGO 文件。用一個(gè)模板文件template.lng把要掃描的參數(shù)寫成占位值再用 Python 循環(huán)替換調(diào)lingo.exe一次跑一個(gè)版本。import subprocess import re rates [0.05, 0.10, 0.15] for i, rate in enumerate(rates): with open(template.lng, r, encodingutf-8) as f: text f.read() text text.replace(DISCOUNT 0;, fDISCOUNT {rate};) with open(frun_{i}.lng, w, encodingutf-8) as f: f.write(text) subprocess.run([lingo.exe, frun_{i}.lng], capture_outputTrue, textTrue, timeout300) with open(lingsol.txt, r, encodingutf-8) as f: out f.read() obj re.search(rObjective value:\s([\d.]), out).group(1) print(frate{rate}, objective{obj})這套流程我在某定價(jià)模型上跑過 20 組參數(shù)全部自動(dòng)落盤不用守著求解器點(diǎn)按鈕。有一點(diǎn)要提醒不同 LINGO 版本命令行輸出的重定向行為不完全一樣如果lingsol.txt沒生成先在 Options 里把 Output Level 調(diào)到 Verbose再確認(rèn)當(dāng)前目錄權(quán)限。把這一套從集合建模、求解器選型、報(bào)告解讀到批量聯(lián)調(diào)走通之后你會(huì)發(fā)現(xiàn)大部分卡殼不是數(shù)學(xué)不會(huì)而是建模習(xí)慣沒改過來。那次我把 Gap 設(shè)成 0 硬跑十萬變量 MIP四個(gè)小時(shí)沒出結(jié)果從那以后我每建一個(gè)新模型都強(qiáng)制先過三件事變量整數(shù)聲明寫完沒、求解器狀態(tài)是 Global 還是 Local、Range 表里要用的影子價(jià)格在不在允許范圍內(nèi)。這幾個(gè)習(xí)慣救了我很多次希望也能幫到你。本文還有配套的精品資源點(diǎn)擊獲取 LINGO建模實(shí)戰(zhàn)：從集合語法到求解器選型與靈敏度分析
簡(jiǎn)介L(zhǎng)INGO線性與非線性全局優(yōu)化器是一款功能強(qiáng)大的數(shù)學(xué)建模優(yōu)化軟件本資源提供其破解安裝版本面向工程、經(jīng)濟(jì)、金融、運(yùn)籌學(xué)等領(lǐng)域的研究者、工程師及建模競(jìng)賽參與者用于求解線性、非線性與整數(shù)規(guī)劃問題。軟件內(nèi)置專門的建模語言用戶可通過簡(jiǎn)潔語句定義變量、目標(biāo)函數(shù)和約束快速搭建大規(guī)模優(yōu)化模型求解器性能強(qiáng)勁支持混合整數(shù)規(guī)劃和結(jié)果敏感性分析還提供數(shù)據(jù)接口與報(bào)告生成功能便于多場(chǎng)景建模與對(duì)比。壓縮包為rar格式大小約15.9MB輕量易獲取目前已有650人學(xué)習(xí)/下載適合從入門到進(jìn)階的不同層級(jí)用戶。借助此資源可直接體驗(yàn)LINGO 11版本的增強(qiáng)求解能力與界面優(yōu)化結(jié)合官方文檔和實(shí)例在供應(yīng)鏈調(diào)度、投資組合優(yōu)化、能源系統(tǒng)設(shè)計(jì)等實(shí)際場(chǎng)景中快速構(gòu)建最優(yōu)解方案提升建模與決策效率。1. LINGO 不是計(jì)算器它到底在解什么LINGO 這個(gè)工具表面看是把數(shù)學(xué)公式敲進(jìn)去、點(diǎn) Solve 就出結(jié)果的“計(jì)算器”但真正拉開建模差距的是它對(duì)集合、數(shù)據(jù)段和求解器分組這三個(gè)設(shè)計(jì)點(diǎn)的處理。同一道排產(chǎn)題有人用幾十行腳本勉強(qiáng)跑通有人用幾行集合寫法在秒級(jí)拿到最優(yōu)解差別不在數(shù)學(xué)功底而在建模方式。這篇文章不打算按教科書順序講語法而是從建模語法、求解器選型、報(bào)告解讀到批量聯(lián)調(diào)完整走一遍我在真實(shí)項(xiàng)目里驗(yàn)證過的流程。適合手頭有優(yōu)化問題的從業(yè)者——供應(yīng)鏈調(diào)度、生產(chǎn)計(jì)劃、定價(jià)決策都算——也適合想快速驗(yàn)證算法思路的同學(xué)。建議開著 LINGO 邊看邊跑大部分示例改兩行數(shù)據(jù)就能直接復(fù)現(xiàn)。2. 集合與數(shù)據(jù)段把模型從紙面搬進(jìn) LINGO 的兩種寫法2.1 原始集與派生集先把“對(duì)誰建?！睂懬宄﨤INGO 和直接寫 LP/MPS 文件的習(xí)慣很不一樣。從紙面公式搬到 LINGO第一件事不是寫目標(biāo)函數(shù)而是定義集合。集合相當(dāng)于告訴求解器我要對(duì)哪些對(duì)象建模型每個(gè)對(duì)象帶哪些屬性。以最經(jīng)典的運(yùn)輸問題為例需求端是 4 個(gè)客戶供給端是 3 個(gè)倉庫每個(gè)倉庫有供應(yīng)量每個(gè)客戶有需求量每條倉庫到客戶的路徑有單位運(yùn)輸成本和運(yùn)輸量。MODEL: SETS: WAREHOUSE /W1, W2, W3/: SUPPLY; CUSTOMER /C1, C2, C3, C4/: DEMAND; ROUTE(WAREHOUSE, CUSTOMER): COST, VOLUME; ENDSETS DATA: SUPPLY 30 40 50; DEMAND 20 25 30 45; COST 2 4 5 3 3 5 4 6 4 2 3 5; ENDDATA MIN SUM(ROUTE(I, J): COST(I, J) * VOLUME(I, J)); FOR(CUSTOMER(J): SUM(WAREHOUSE(I): VOLUME(I, J)) DEMAND(J)); FOR(WAREHOUSE(I): SUM(CUSTOMER(J): VOLUME(I, J)) SUPPLY(I)); FOR(ROUTE(I, J): VOLUME(I, J) 0); END這段代碼里的三類寫法建議記牢。第一原始集定義用名字 /成員/ : 屬性成員之間用逗號(hào)或空格隔開屬性跟在冒號(hào)后面一個(gè)屬性就是一個(gè)和集合同長(zhǎng)度的參數(shù)或變量數(shù)組。第二派生集ROUTE(WAREHOUSE, CUSTOMER)表示兩個(gè)原始集的笛卡爾積展開后共 12 個(gè)成員COST和VOLUME都掛在這個(gè)派生集上分別讀入 12 個(gè)成本數(shù)據(jù)和生成 12 個(gè)決策變量。第三FOR和SUM是 LINGO 的循環(huán)與求和函數(shù)FOR(CUSTOMER(J): ...)表示對(duì)每個(gè)客戶單獨(dú)生成一條等式約束SUM(WAREHOUSE(I): ...)表示對(duì)倉庫下標(biāo)求和。我一般會(huì)先在草稿紙上畫一張“對(duì)象關(guān)系圖”把集合畫成方框、屬性畫成方框里的行、跨集合的連線就是派生集。這樣做直接的好處是寫FOR時(shí)不容易搞混下標(biāo)。比如上面第二個(gè)FOR里外層循環(huán)下標(biāo)是I倉庫內(nèi)部求和下標(biāo)是J客戶意思是每個(gè)倉庫發(fā)出去的總量等于它的供應(yīng)量如果把兩者位置對(duì)調(diào)模型含義就從“運(yùn)出量守恒”變成“運(yùn)入量等于別的東西”LINGO 不會(huì)報(bào)警告只有結(jié)果跟業(yè)務(wù)對(duì)不上時(shí)你才會(huì)發(fā)現(xiàn)。2.2 數(shù)據(jù)段與外部文件把參數(shù)和模型結(jié)構(gòu)分開上一節(jié)的運(yùn)輸問題里數(shù)據(jù)全寫在DATA/ENDDATA段中這是我最推薦的起始方式模型結(jié)構(gòu)一眼能看完改數(shù)據(jù)也方便。但數(shù)據(jù)量超過幾十個(gè)、或者同一套模型要套用到多個(gè)區(qū)域時(shí)就應(yīng)該把參數(shù)挪到外部文件。LINGO 的FILE函數(shù)就是干這個(gè)的。DATA: SUPPLY FILE(sup_data.txt); DEMAND FILE(dem_data.txt); COST FILE(cost_data.txt); ENDDATAsup_data.txt里只需要按集合成員順序?qū)憯?shù)值例如30 40 50文件里可以有換行和空格數(shù)值之間不需要逗號(hào)。這樣模型文件里只剩集合和約束結(jié)構(gòu)數(shù)據(jù)歸數(shù)據(jù)。我做過的一次模擬項(xiàng)目 X就是同一套庫存補(bǔ)貨模型把華東、華南、華北三組數(shù)據(jù)分別存成三個(gè)文本跑三次模型一行沒改。把數(shù)據(jù)外置還有一個(gè)隱藏優(yōu)點(diǎn)做參數(shù)掃描時(shí)你只需要用腳本替換文本里的某個(gè)值然后重新求解。后面第 6 章我會(huì)給完整的批量做法。如果你想讓結(jié)果導(dǎo)出成可讀文件DATA段里也可以寫TEXT(result.txt)把某個(gè)集合屬性的值按成員順序?qū)懙轿谋疚募?。這個(gè)函數(shù)我在前面的選址問題里常用至少比從密密麻麻的 Solution Report 里手工復(fù)制要穩(wěn)。2.3 從紙面公式到 LINGO 表達(dá)的三個(gè)習(xí)慣最后補(bǔ)三個(gè)我踩過坑練出來的習(xí)慣。第一目標(biāo)函數(shù)寫在最前面MIN或MAX開頭接下來所有等號(hào)、不等式、BND都是約束。第二約束盡量加行標(biāo)簽比如FOR(CUSTOMER(J): [DEMAND_BALANCE] SUM(...) DEMAND(J));加標(biāo)簽之后報(bào)告里的 Row 編號(hào)就會(huì)帶著名字排查結(jié)果時(shí)能少一半時(shí)間去對(duì)行號(hào)。第三決策變量和參數(shù)不要混在一個(gè)集合屬性里參數(shù)在DATA段賦值變量只聲明不賦值如果非要給變量初值可以在DATA段賦但別把初值當(dāng)業(yè)務(wù)參數(shù)用否則下次換數(shù)據(jù)時(shí)容易漏改。3. 求解器選型LP、MIP、NLP 跑不對(duì)時(shí)先查這三處3.1 線性模型與整數(shù)模型變量類型聲明和 Gap 設(shè)置LINGO 會(huì)根據(jù)模型自動(dòng)分派求解器純線性連續(xù)問題是 LP帶整數(shù)變量是 MIP目標(biāo)或約束里有非線性表達(dá)式就是 NLP。自動(dòng)分派省事但很多“跑出來不對(duì)”的翻車現(xiàn)場(chǎng)恰恰是因?yàn)槟Ｐ捅蛔詣?dòng)分派到了錯(cuò)誤的求解器。最典型的是整數(shù)變量沒聲明模型被當(dāng)成 LP 解得到小數(shù)方案。設(shè)施選址問題能同時(shí)展示 LP 到 MIP 的差別。MODEL: SETS: PLANT /P1, P2, P3/: FIXED, CAPACITY, OPEN; CUSTOMER /C1, C2, C3, C4/: DEMAND; ROUTE(PLANT, CUSTOMER): TRANSCOST, SHIP; ENDSETS DATA: FIXED 1000 800 1200; CAPACITY 500 400 600; DEMAND 150 180 120 200; TRANSCOST 5 4 6 3 7 3 5 6 4 5 4 7; ENDDATA MIN SUM(PLANT(I): FIXED(I) * OPEN(I)) SUM(ROUTE(I, J): TRANSCOST(I, J) * SHIP(I, J)); FOR(CUSTOMER(J): SUM(PLANT(I): SHIP(I, J)) DEMAND(J)); FOR(PLANT(I): SUM(CUSTOMER(J): SHIP(I, J)) CAPACITY(I) * OPEN(I)); FOR(ROUTE(I, J): SHIP(I, J) 0); FOR(PLANT(I): BIN(OPEN(I))); END關(guān)鍵在最后兩行。BIN(OPEN(I))把OPEN聲明成 0/1 變量模型才是 MIP倒數(shù)第二條容量約束 CAPACITY(I) * OPEN(I)的含義是OPEN為 1 時(shí)該廠產(chǎn)能正常生效為 0 時(shí)上限變成 0運(yùn)輸量被迫為 0。這個(gè)“乘法開關(guān)”模式幾乎適用于所有帶固定成本的選擇問題。MIP 跑得慢或者太久不出最優(yōu)解時(shí)先別急著加機(jī)器去 Options 里的 MIP Solver 參數(shù)看 Dual Gap。Gap 表示當(dāng)前最好的整數(shù)可行解與松弛上界之間的相對(duì)差距設(shè)成 0 表示一定要證明最優(yōu)才停代價(jià)是分支定界可能跑幾個(gè)小時(shí)設(shè)成 1% 通常能在十幾分鐘內(nèi)給出工程上可接受的解。我在做某跨平臺(tái)系統(tǒng)的排產(chǎn)模塊時(shí)就把 Gap 設(shè)成 0.5%業(yè)務(wù)方完全接受求解時(shí)間從 40 分鐘壓到 3 分鐘。要注意如果模型本身只有幾千個(gè)變量Gap 設(shè)大點(diǎn)無所謂一旦變量上十萬Gap 開 1% 和開 5% 的耗時(shí)差能到 10 倍以上量力設(shè)。3.2 非線性模型與全局求解器默認(rèn)只是局部收斂非線性模型的默認(rèn)行為容易讓人誤判。LINGO 碰到非線性目標(biāo)或約束時(shí)默認(rèn)走 NLP 求解器本質(zhì)是局部搜索從一個(gè)初始點(diǎn)出發(fā)沿可行方向逐步改進(jìn)找到一個(gè)“站在附近看確實(shí)最優(yōu)”的解就停。如果模型是凸的局部即全局但現(xiàn)實(shí)中定價(jià)、選址、排程模型很少有凸的于是局部解就變成一個(gè)大坑。一個(gè)能快速?gòu)?fù)現(xiàn)的示例是圓柱儲(chǔ)罐用料最小化給定容積 100 立方米求半徑 R 和高度 H 使表面積最小。MODEL: ! 目標(biāo)最小化圓柱儲(chǔ)罐的用料表面積 V 100; ! 容量固定 100 立方米 MIN 2 * 3.14159 * R * R 2 * 3.14159 * R * H; BND(1, R, 10); BND(1, H, 50); 3.14159 * R * R * H V; END這個(gè)模型手算就能驗(yàn)證把 H 換成 100 / (3.14159 * R^2)目標(biāo)函數(shù)變成一元函數(shù)最優(yōu)在 R 約 2.5、H 約 5.0材料面積大概 118 平方米。默認(rèn) NLP 從初值 R1、H20 也能找到這個(gè)點(diǎn)因?yàn)榭尚杏蚯嫔线@個(gè)點(diǎn)是唯一谷底。但如果把目標(biāo)函數(shù)換成帶SIN或EXP的生產(chǎn)計(jì)劃模型谷底可能就不止一個(gè)。這時(shí)候要在 LINGO 菜單里打開 Options→Global Solver勾選啟用全局求解器它會(huì)在后臺(tái)做類似分支定界的過程反復(fù)逼近并證明全局最優(yōu)。代價(jià)是時(shí)間變量一多、函數(shù)一復(fù)雜跑幾十分鐘甚至幾天都正常。我一般先把局部解跑出來看結(jié)果是否在業(yè)務(wù)合理區(qū)間再用全局求解器驗(yàn)證一次而不是一開始就全局硬算。3.3 三個(gè)參數(shù)與數(shù)值尺度TIMLIM、Scale Model 和線性化選項(xiàng)求解器參數(shù)里我最常動(dòng)的是三個(gè)。第一個(gè)是時(shí)間上限 TIMLIM單位秒設(shè)到 300 到 600 很常見它防止 MIP 或全局求解器無限跑下去。第二個(gè)是迭代上限 ITERLIM只對(duì)單純形法有意義內(nèi)點(diǎn)法和分支定界主要看時(shí)間。第三個(gè)是 Scale Model 選項(xiàng)遇到目標(biāo)上百萬、約束系數(shù)只有 0.001 這類數(shù)值尺度失衡時(shí)勾上它經(jīng)常會(huì)讓 Infeasibilities 從 0.01 突然變成 0.000000。原理很簡(jiǎn)單矩陣?yán)锿瑫r(shí)存在 1e6 和 1e-4 級(jí)的數(shù)求解器的內(nèi)積、求逆精度會(huì)吃緊縮放后條件數(shù)改善結(jié)果立刻穩(wěn)。另外 LINGO 有自動(dòng)線性化能力小規(guī)模非線性表達(dá)式會(huì)在預(yù)求解階段嘗試線性化但一旦變量超過幾百個(gè)就別指望它全自動(dòng)人工引入輔助變量改寫反而更快這部分在第 5 章會(huì)給改寫例子。4. 報(bào)告與靈敏度分析從 Solution Report 里讀出真實(shí)邊界4.1 先看懂求解狀態(tài)和三個(gè)關(guān)鍵字段每次求解完彈出的 Solution Report大多數(shù)人只掃一眼目標(biāo)值就關(guān)掉。我建議至少看四個(gè)地方。首先是報(bào)告第一行Global optimal solution found說明已經(jīng)證明全局最優(yōu)Local optimum found說明只是局部解Feasible solution found更弱只保證找到可行解。這三者差一個(gè)級(jí)別業(yè)務(wù)決策的置信度就差一個(gè)檔次。其次看Infeasibilities它不為 0 說明模型找不到嚴(yán)格滿足所有約束的點(diǎn)這時(shí)目標(biāo)值再好看都不能用。下面是一個(gè)簡(jiǎn)化后的 LP 小案例報(bào)告用來解釋字段含義不必對(duì)應(yīng)前面的某個(gè)具體模型。Global optimal solution found. Objective value: 870.0000 Objective bound: 870.0000 Infeasibilities: 0.000000 Total solver iterations: 5 Variable Value Reduced Cost X1 40.00000 0.000000 X2 20.00000 5.000000 Row Slack or Surplus Dual Price 1 870.0000 1.000000 2 0.000000 7.000000 3 10.00000 0.000000變量區(qū)里的Reduced Cost對(duì)基變量通常為 0對(duì)非基變量表示如果把這個(gè)變量強(qiáng)行加進(jìn)最優(yōu)解一個(gè)單位目標(biāo)值會(huì)惡化多少。Row 區(qū)里的Slack or Surplus是約束的剩余量等于 0 說明該約束緊貼邊界是真正的限制因素大于 0 說明約束還有余量暫時(shí)不影響目標(biāo)。Dual Price就是影子價(jià)格下面單獨(dú)說。Row 后面的編號(hào)對(duì)應(yīng)模型里的約束順序所以我在第 2.3 節(jié)強(qiáng)調(diào)要加行標(biāo)簽就是為了讓這里的編號(hào)能一眼對(duì)上約束名。我在實(shí)戰(zhàn)里的習(xí)慣是改動(dòng)模型后先看求解狀態(tài)那行字再看 Infeasibilities最后才看目標(biāo)值。狀態(tài)不對(duì)時(shí)數(shù)值再漂亮也只是黑匣子輸出。4.2 靈敏度 Range 表影子價(jià)格能走多遠(yuǎn)默認(rèn)報(bào)告不會(huì)給完整的靈敏度 Range 表。要在求解前把 Options 里 General Solver 的 Dual Computations 設(shè)成 Prices Ranges重新求解后才會(huì)出現(xiàn)兩組區(qū)間目標(biāo)系數(shù)變化范圍和右側(cè) RHS 變化范圍。目標(biāo)系數(shù)范圍表長(zhǎng)這樣。變量當(dāng)前系數(shù)允許增加允許減少X150.0010.0020.00X240.008.006.00RHS 范圍表單獨(dú)一張。約束當(dāng)前 RHS允許增加允許減少約束2100.005.008.00約束380.0012.0020.00這兩張表回答的問題是我不重新求解的情況下業(yè)務(wù)參數(shù)在多大范圍內(nèi)波動(dòng)當(dāng)前這組最優(yōu)決策仍然是“結(jié)構(gòu)最優(yōu)”。注意“結(jié)構(gòu)最優(yōu)”不等于“變量值不變”它只表示最優(yōu)基的變量集合不變。一旦某個(gè)系數(shù)的實(shí)際增幅超出允許增加最優(yōu)基會(huì)變影子價(jià)格會(huì)跳到另一個(gè)值這時(shí)候再去引用原 Dual Price 就作廢了。影子價(jià)格業(yè)務(wù)化是最容易翻車的地方。某一次我給一個(gè)倉儲(chǔ)項(xiàng)目做資源評(píng)估報(bào)告里約束 2 的 Dual Price 是 7意思是該約束右側(cè)再給 1 個(gè)單位目標(biāo)值能改善 7。業(yè)務(wù)方很興奮說多放 10 個(gè)貨架就多賺 70 萬。我拿 Range 表一看這個(gè)約束的 RHS 允許增加只有 5第 6 個(gè)貨架開始 Shadow Price 已經(jīng)不是 7 了。所以我把影子價(jià)格當(dāng)“邊際價(jià)值信號(hào)”不當(dāng)“線性外推的乘法口訣表”并且只對(duì)它生效的范圍負(fù)責(zé)。4.3 用手工可解案例做回歸報(bào)告說最優(yōu)不等于邏輯對(duì)還有一層驗(yàn)證是很多人跳過的報(bào)告說最優(yōu)不等于模型建得沒問題。我習(xí)慣在正式跑大規(guī)模數(shù)據(jù)前先構(gòu)造一個(gè)規(guī)模小到能手工算的案例塞進(jìn)同一套模型確認(rèn)結(jié)果和手算一致。比如運(yùn)輸問題我做過 2 個(gè)倉庫送 2 個(gè)客戶的小案例手算最優(yōu)值再讓模型跑對(duì)上了才換真實(shí)數(shù)據(jù)。這件事花 10 分鐘能擋住一大半“集合寫錯(cuò)但報(bào)告照出最優(yōu)解”的隱性 bug。5. 避坑實(shí)錄五個(gè)反復(fù)出現(xiàn)的建模翻車場(chǎng)景與處置5.1 整數(shù)聲明缺失解出 13.7 臺(tái)設(shè)備業(yè)務(wù)當(dāng)場(chǎng)不認(rèn)現(xiàn)象生產(chǎn)計(jì)劃模型跑完報(bào)告里生產(chǎn)量是 13.7 臺(tái)、開啟 0.4 條產(chǎn)線業(yè)務(wù)方直接拒收結(jié)果。原因LINGO 變量默認(rèn)是連續(xù)變量GIN、BIN沒寫模型被當(dāng)成 LP 解了。解決凡是決策量必須是整數(shù)的聲明時(shí)給集合屬性補(bǔ)一行GIN(PRODUCE(P))或BIN(OPEN(I))。跑完先看報(bào)告底部的 Model Statistics核對(duì)整數(shù)變量個(gè)數(shù)。我后來養(yǎng)成的習(xí)慣是需求評(píng)審階段就把“哪些量必須是整數(shù)”寫進(jìn)建模文檔避免寫代碼時(shí)漏掉。5.2 集合維度錯(cuò)配Error 72 的鍋十有八九在這現(xiàn)象加載模型報(bào) Runtime Error錯(cuò)誤碼 72指向某一行FOR或SUM。原因下標(biāo)集合寫錯(cuò)了維度。比如在FOR(CUSTOMER(J): ...)內(nèi)部卻用了SUM(WAREHOUSE(J): SHIP(I, J))把客戶集合的成員J當(dāng)成了倉庫集合的成員用。解決先定位報(bào)錯(cuò)行檢查每個(gè)循環(huán)變量所屬的集合。FOR(A(I): ...)里可以安全使用的下標(biāo)只有I和它所在約束里通過SUM跨集合引用的其他集合下標(biāo)。一個(gè)笨但有效的辦法把模型里所有FOR和SUM的下標(biāo)變量名列出來逐個(gè)對(duì)照 SETS 定義。另外派生集下標(biāo)的順序要和定義時(shí)一致ROUTE(WAREHOUSE, CUSTOMER)里第一維是倉庫第二維是客戶寫ROUTE(J, I)卻把 I 當(dāng)倉庫、J 當(dāng)客戶就會(huì)維度錯(cuò)配。5.3 IF 和 SMAX 把問題拆成非光滑懸崖現(xiàn)象模型能解但換一組合理初值結(jié)果跳到一個(gè)完全不同的解業(yè)務(wù)怎么調(diào)參數(shù)都達(dá)不到報(bào)告里的“最優(yōu)”。原因目標(biāo)或約束里用了IF、SMAX、ABS這類函數(shù)可行域被切成不連續(xù)的多塊非凸不可導(dǎo)默認(rèn) NLP 局部搜索只能爬到其中一個(gè)谷底。解決把非光滑項(xiàng)改寫成線性不等式。比如約束里出現(xiàn)SMAX(0, PRICE(I) - BID(I))可以引入輔助變量OVER(I)替換成下面兩條約束。OVER(I) PRICE(I) - BID(I); OVER(I) 0;原來用SMAX的地方直接換成OVER(I)。數(shù)學(xué)上等價(jià)數(shù)值上卻從“不可導(dǎo)懸崖”變成“兩條線性邊界”全局求解器才真正能用。類似地IF(COND, A, B)可以分解成一個(gè) 0/1 變量加兩個(gè)大 M 不等式。需要提醒的是大 M 不要拍腦袋給 1e6取“業(yè)務(wù)上該量級(jí)最大可能值的 3 倍”就夠太大反而把數(shù)值條件數(shù)拖壞。5.4 數(shù)值尺度懸殊手算滿足約束報(bào)告卻報(bào) Infeasibilities現(xiàn)象目標(biāo)函數(shù)動(dòng)輒上百萬約束系數(shù)小到 0.0001報(bào)告里 Infeasibilities 顯示 0.01可人工逐個(gè)約束驗(yàn)算發(fā)現(xiàn)都滿足。原因模型矩陣的條件數(shù)太差求解器在數(shù)值容差內(nèi)判斷違反約束產(chǎn)生了“數(shù)值偽不可行”。這在小數(shù)單位和萬元單位混用的成本核算模型里很常見。解決先統(tǒng)一業(yè)務(wù)單位把元換成萬元或把千克換成噸讓目標(biāo)值和系數(shù)都落在 1e-2 到 1e4 之間再去 Options 里勾選 Scale Model。我在某物流成本模型里遇到過 Infeasibilities0.03把單箱成本從“元”改成“萬元”后直接歸零一秒沒多算。如果改完還不行重點(diǎn)查是否有同時(shí)除以兩個(gè)量級(jí)相差很大的分母的表達(dá)式。5.5 局部最優(yōu)當(dāng)全局最優(yōu)同一個(gè)模型換初值換結(jié)果現(xiàn)象同一個(gè) NLP 模型同事用另一套求解器算出來的目標(biāo)值明顯比我的好換初始點(diǎn)再跑結(jié)果又變一個(gè)樣。原因默認(rèn) NLP 是局部搜索非凸問題從不同起點(diǎn)收斂到不同局部解根本沒有“唯一答案”。解決重要決策用全局求解器跑一遍同時(shí)設(shè)時(shí)間上限比如 600 秒看它找到的目標(biāo)值和默認(rèn) NLP 是否一致。如果全局求解時(shí)間完全不可接受退一步做多起點(diǎn)掃描人工挑 5 到 10 個(gè)有業(yè)務(wù)含義的初始點(diǎn)跑局部 NLP記錄每個(gè)起點(diǎn)對(duì)應(yīng)的目標(biāo)值和決策量取業(yè)務(wù)上最合理的那組并在文檔里注明這是“驗(yàn)證過的局部解”不是全局最優(yōu)。從那以后我對(duì)“最優(yōu)”這個(gè)詞都加限定語是全局證明還是局部搜索加經(jīng)驗(yàn)判斷。5.6 排查速查表癥狀最先查的語法/參數(shù)本章對(duì)應(yīng)場(chǎng)景解出小數(shù)GIN / BIN 聲明5.1Runtime Error 72FOR/SUM 的下標(biāo)集合匹配5.2換初值結(jié)果大變IF/SMAX/ABS 非線性改寫5.3Infeasibilities 不歸零Scale Model、單位歸一化5.4報(bào)告只有 Local optimum開 Global Solver5.56. 把 LINGO 嵌進(jìn)工作流Excel 聯(lián)調(diào)與批量參數(shù)掃描6.1 用 OLE 讓 Excel 當(dāng)數(shù)據(jù)界面Excel 里先給區(qū)域命名LINGO 的OLE就能按名字讀寫。數(shù)據(jù)段讀入?yún)?shù)CALC 段寫回結(jié)果業(yè)務(wù)方改 Excel 就行不碰模型文件。DATA: SUPPLY OLE(C:\supply_plan.xlsx, supply); DEMAND OLE(C:\supply_plan.xlsx, demand); COST OLE(C:\supply_plan.xlsx, cost); ENDDATA CALC: OLE(C:\supply_plan.xlsx, ship_out) SHIP; ENDCALC注意路徑里別用中文名Excel 文件要事先打開或在名稱管理器里確認(rèn)區(qū)域存在。第一次聯(lián)調(diào)失敗九成是 Excel 里的命名區(qū)域尺寸和集合成員數(shù)不一致。提示Excel 命名區(qū)域的尺寸必須和集合成員數(shù)一致否則 OLE 讀入會(huì)在運(yùn)行時(shí)報(bào)告維度不匹配先把區(qū)域重新框選再命名一次。6.2 參數(shù)掃描改文件、跑求解、讀結(jié)果需要做敏感性掃描時(shí)我不喜歡手動(dòng)改 LINGO 文件。用一個(gè)模板文件template.lng把要掃描的參數(shù)寫成占位值再用 Python 循環(huán)替換調(diào)lingo.exe一次跑一個(gè)版本。import subprocess import re rates [0.05, 0.10, 0.15] for i, rate in enumerate(rates): with open(template.lng, r, encodingutf-8) as f: text f.read() text text.replace(DISCOUNT 0;, fDISCOUNT {rate};) with open(frun_{i}.lng, w, encodingutf-8) as f: f.write(text) subprocess.run([lingo.exe, frun_{i}.lng], capture_outputTrue, textTrue, timeout300) with open(lingsol.txt, r, encodingutf-8) as f: out f.read() obj re.search(rObjective value:\s([\d.]), out).group(1) print(frate{rate}, objective{obj})這套流程我在某定價(jià)模型上跑過 20 組參數(shù)全部自動(dòng)落盤不用守著求解器點(diǎn)按鈕。有一點(diǎn)要提醒不同 LINGO 版本命令行輸出的重定向行為不完全一樣如果lingsol.txt沒生成先在 Options 里把 Output Level 調(diào)到 Verbose再確認(rèn)當(dāng)前目錄權(quán)限。把這一套從集合建模、求解器選型、報(bào)告解讀到批量聯(lián)調(diào)走通之后你會(huì)發(fā)現(xiàn)大部分卡殼不是數(shù)學(xué)不會(huì)而是建模習(xí)慣沒改過來。那次我把 Gap 設(shè)成 0 硬跑十萬變量 MIP四個(gè)小時(shí)沒出結(jié)果從那以后我每建一個(gè)新模型都強(qiáng)制先過三件事變量整數(shù)聲明寫完沒、求解器狀態(tài)是 Global 還是 Local、Range 表里要用的影子價(jià)格在不在允許范圍內(nèi)。這幾個(gè)習(xí)慣救了我很多次希望也能幫到你。本文還有配套的精品資源點(diǎn)擊獲取