簡介線性規(guī)劃單純形算法的C實現(xiàn)資料包面向運籌學(xué)課程學(xué)習(xí)者、算法愛好者及需要求解線性規(guī)劃問題的開發(fā)者演示如何將標(biāo)準(zhǔn)線性規(guī)劃模型轉(zhuǎn)化為單純形表并迭代求最優(yōu)解。壓縮包共17個文件、大小623KB包含兩個C源文件及頭文件、Visual C 6.0工程文件dsp/dsw/ncb/opt/plg、編譯調(diào)試產(chǎn)物obj/pdb/ilk/idb/exe以及輸入輸出示例文本可直接閱讀源碼或運行程序驗證。已有1050人學(xué)習(xí)下載適合用于課程設(shè)計、算法實現(xiàn)參考或自學(xué)單純形法原理。從標(biāo)準(zhǔn)輸入文件讀取數(shù)據(jù)經(jīng)歷建立模型、構(gòu)建初始單純形表、檢驗數(shù)判斷、基變量替換等流程讀者可對照源碼理解每一步驟尤其是單純形表更新與檢驗數(shù)判斷的具體實現(xiàn)并通過可執(zhí)行文件快速查看求解結(jié)果附帶的工程文件也便于二次修改與調(diào)試是一份緊湊實用的算法學(xué)習(xí)樣例。 很多人在學(xué)算法時會把線性規(guī)劃當(dāng)成一門純理論課覺得它離業(yè)務(wù)很遠(yuǎn)。我過去也這么想直到第一次做排產(chǎn)優(yōu)化時被現(xiàn)實教育了一頓——生產(chǎn)計劃、物流派車、人力排班、庫存?zhèn)湄涍@些業(yè)務(wù)問題吵來吵去其實都能落成同一套數(shù)學(xué)模型而處理這套模型的算法就是線性規(guī)劃。換句話說線性規(guī)劃不是考試用完了就扔的公式它是那種“你越早會用越早受益”的決策工具。這篇文章我會從問題定義、數(shù)學(xué)直覺、建模方法、Python工具鏈、完整案例到常見坑位全部過一遍特別適合后端、算法、數(shù)據(jù)崗位的同學(xué)以及所有需要做資源分配決策但不想每次靠拍腦袋的人。讀完你至少能做到拿到一個業(yè)務(wù)需求能判斷能不能用線性規(guī)劃能獨立建出模型能在求解器報錯時知道問題出在哪兒。1. 線性規(guī)劃到底在解決什么問題不止是“求最優(yōu)解”這么簡單1.1 線性規(guī)劃問題的三要素任何一個線性規(guī)劃問題拆到最底層只有三樣?xùn)|西決策變量、目標(biāo)函數(shù)、約束條件。決策變量是你要拍板的值比如“A產(chǎn)品生產(chǎn)多少件”“給B城市發(fā)幾車貨”目標(biāo)函數(shù)是一個關(guān)于決策變量的線性表達式用來衡量方案好壞通常寫成最大化利潤或最小化成本約束條件則是一組線性不等式或等式描述資源上限、需求下限、工序先后等現(xiàn)實限制。這里“線性”兩個字很關(guān)鍵。目標(biāo)函數(shù)和約束條件里的變量只能是一次方不能出現(xiàn) x2、sin(x)、x*y 這類非線性項。原因在于一旦非線性問題的幾何性質(zhì)會完全改變后面要講的單純形法也不再適用?？梢源直┑乩斫獬删€性保證了“效果可預(yù)期、結(jié)果可計算”這是它能在工業(yè)界大規(guī)模落地的根本原因。1.2 什么樣的問題適合用線性規(guī)劃我判斷一個業(yè)務(wù)問題能不能用線性規(guī)劃基本就看四點決策結(jié)果能不能用一組實數(shù)或整數(shù)變量描述優(yōu)化目標(biāo)能不能寫成決策變量的加權(quán)和限制條件能不能寫成線性不等式你要的是不是全局最優(yōu)而不是“經(jīng)驗上差不多”。如果這四條的答案都是“是”那大概率可以轉(zhuǎn)成線性規(guī)劃。制造業(yè)排產(chǎn)、物流配送路徑選擇、倉儲補貨、電力調(diào)度、投資組合配置這些場景表面差別很大模型結(jié)構(gòu)卻很相似差別只在變量和約束的數(shù)量與含義。1.3 它和普通算法題最大的差異很多人學(xué)數(shù)據(jù)結(jié)構(gòu)時習(xí)慣了“給定輸入求輸出”的思路排序、搜索、動態(tài)規(guī)劃都是這個模式核心是設(shè)計計算過程。線性規(guī)劃不一樣它的輸入是一套規(guī)則和限制核心是“怎么從無數(shù)種可行方案里挑一個最優(yōu)決策”。這個差異決定了思維方式必須從“寫計算邏輯”切換成“寫業(yè)務(wù)約束”。我見過不少科班出身的人拿到業(yè)務(wù)第一步就想著要不要用循環(huán)、遞歸結(jié)果繞了一大圈其實用建模語言把約束寫清楚求解器幾秒鐘就出結(jié)果了。順帶辟個謠線性規(guī)劃和線性回歸是兩碼事。前者是數(shù)學(xué)規(guī)劃里的優(yōu)化模型后者是統(tǒng)計學(xué)里擬合數(shù)據(jù)的工具名字長得像解決的是完全不同的兩類問題。2. 數(shù)學(xué)直覺與單純形法最優(yōu)解為什么總是“跑在邊界上”2.1 二維情形的幾何直覺先看只有兩個決策變量的情況這是最好理解的。每個線性約束對應(yīng)平面上的一條直線直線把平面切成兩個半平面所有約束相交出來的區(qū)域是一個凸多邊形這個多邊形就是可行域。目標(biāo)函數(shù) z 3x 4y 在固定 z 值時是一條直線我們做的事相當(dāng)于把這條直線沿著利潤增大的方向平移直到它剛好還在可行域上碰到某個極限位置。這個極限位置永遠(yuǎn)是多邊形的頂點不會是邊中間更不會在內(nèi)部。你可以想象一個略微傾斜的桌面桌上放著一個多邊形托盤托盤里有一顆珠子桌面整體是平坦傾斜的珠子最后停的位置一定是托盤邊緣的某個拐角不會懸在中間。這就是線性規(guī)劃最優(yōu)解總是落在頂點上的直覺來源。2.2 單純形法從一個頂點走到另一個頂點單純形法正是利用了上面這個性質(zhì)。它的核心思想是從可行域的一個頂點出發(fā)沿著某條棱邊走到相鄰頂點如果這個頂點的目標(biāo)函數(shù)值更好就繼續(xù)換直到找不到更好的相鄰頂點為止。這個過程和爬山很像。從山腳出發(fā)沿著山脊往上爬每到一個埡口就看看四周有沒有更高的點有就繼續(xù)走沒有就說明到山頂了。單純形法的精妙之處在于它不需要遍歷所有頂點——實際問題里頂點數(shù)量可能多到爆炸但只要沿著能讓目標(biāo)改善的方向走通常幾十步甚至幾步就能收斂。2.3 內(nèi)點法和整數(shù)規(guī)劃的一筆帶過單純形法雖然經(jīng)典但遇到特別大規(guī)模的問題時可能需要頻繁變換基變量性能會受影響。于是就有了內(nèi)點法不沿著邊界走而是直接從可行域內(nèi)部向最優(yōu)頂點逼近。現(xiàn)代求解器比如 HiGHS、Gurobi 通常都會內(nèi)置多種算法根據(jù)問題特征自動切換用戶基本不用關(guān)心底層用的是什么。另外要提一句整數(shù)規(guī)劃。很多實際問題不允許變量取小數(shù)比如“派幾輛車”不可能是3.7輛。這類問題叫整數(shù)規(guī)劃或混合整數(shù)規(guī)劃MILP求解方法是在線性規(guī)劃的基礎(chǔ)上做分支定界。所以想玩轉(zhuǎn)整數(shù)規(guī)劃先把線性規(guī)劃搞明白是必須的否則連門都摸不到。3. 建模的思維方式把現(xiàn)實約束翻譯成數(shù)學(xué)語言3.1 從命令式思維到聲明式思維程序員最大的坎往往不是數(shù)學(xué)而是思維轉(zhuǎn)換。寫算法題時習(xí)慣了命令式思維一步一步告訴機器“先做這個再做那個”線性規(guī)劃要求的是聲明式思維你只負(fù)責(zé)把業(yè)務(wù)規(guī)則翻譯成數(shù)學(xué)式子至于怎么求最優(yōu)解那是求解器的事。打個比方命令式思維是“你親自開車每個路口都要判斷怎么走”聲明式思維是“你告訴司機目的地和不能走的路剩下交給他”。初學(xué)者最常犯的錯誤就是在模型里試圖教求解器“應(yīng)該怎么算”結(jié)果把模型搞得一團糟。3.2 建模五步法我自己建模有一套固定流程每次照著走能省下大量返工時間列出所有決策變量寫清楚每個變量的含義和單位寫出目標(biāo)函數(shù)確認(rèn)是最大化還是最小化逐條列出業(yè)務(wù)限制翻譯成線性不等式或等式補上默認(rèn)約束最常見的是變量非負(fù)以及變量的上下限跑求解器根據(jù)結(jié)果反推模型是否有遺漏或方向錯誤。這套流程看起來簡單但第三步和第四步最容易出錯。漏一條約束模型可能直接無解多一條錯誤約束解出來可能完全不符合業(yè)務(wù)直覺。3.3 一個最簡單的建模示例用一個幾乎不能再小的例子說明。某工廠生產(chǎn) A、B 兩種產(chǎn)品A 每件利潤 7 元需要 2 小時設(shè)備時間B 每件利潤 5 元需要 1 小時設(shè)備時間設(shè)備每天最多 10 小時A 每天最多生產(chǎn) 3 件問最優(yōu)產(chǎn)量是多少。設(shè) x 為 A 產(chǎn)量y 為 B 產(chǎn)量模型就是max Z 7x 5y2x y ≤ 10x ≤ 3x, y ≥ 0手動解一下x 取滿 3剩余設(shè)備 4 小時全部給 y得到 y4總利潤 7×35×441。這個例子小得不能再小但它完整展示了三要素的翻譯過程新手建議從這個粒度開始練手。3.4 新手最容易踩的三個建?？拥谝粋€坑是把非線性邏輯硬寫成線性約束。比如“如果生產(chǎn) A 就必須生產(chǎn) B”這種條件關(guān)系實際應(yīng)該引入 0-1 變量來表達而不是在約束里寫 x*y 這種交叉項。第二個坑是漏掉變量的上下限和非負(fù)約束。少了這些可行域可能變成無界區(qū)域求解器會直接報 “Unbounded”你就得回頭補約束。第三個坑是量綱不統(tǒng)一。我見過一個排產(chǎn)模型利潤按“萬元/噸”算、工時按“小時/噸”算結(jié)果寫約束時把利潤數(shù)值當(dāng)工時系數(shù)填了進去求解器照樣解出“最優(yōu)解”但結(jié)果完全失真。建模完成后一定要回頭檢查每個系數(shù)的單位和含義這是成本最低的一步體檢。4. 從模型到求解器Python工具鏈怎么選、怎么用4.1 主流工具橫向?qū)Ρ热绻阌?Python 做線性規(guī)劃市面上可選工具不少我按實際使用體驗整理了一個對照表工具適用場景上手難度說明Excel 規(guī)劃求解小規(guī)模、一次性分析低內(nèi)置 Solver適合臨時算一下SciPy linprogPython 數(shù)據(jù)生態(tài)的中等規(guī)模 LP中API 簡潔MILP 支持有限PuLP中小規(guī)模 LP 和 MILP 建模低語法接近數(shù)學(xué)表達式推薦入門OR-Tools大規(guī)模 MILP、CP-SAT、排班類問題中高Google 出品功能全面Gurobi / CPLEX企業(yè)級大規(guī)模問題中高商業(yè)許可性能與穩(wěn)定性頂尖4.2 為什么推薦從 PuLP 入手如果只選一個工具入門我會選 PuLP。理由很簡單語法和數(shù)學(xué)表達式幾乎一一對應(yīng)寫出來的代碼可以直接對照公式檢查開源免費pip 一條命令裝完默認(rèn)自帶 CBC 求解器線性規(guī)劃和整數(shù)規(guī)劃都能解覆蓋了大部分業(yè)務(wù)場景。更關(guān)鍵的是PuLP 的模型代碼和后端求解器是解耦的。同一個模型你可以在 CBC、GLPK、HiGHS 甚至 Gurobi 之間切換模型本身不用改。這意味著你不需要在入門階段就綁定某個商業(yè)化產(chǎn)品以后規(guī)模大了再換求解器也來得及。4.3 求解器內(nèi)核到底是干嘛的很多人會把 PuLP 和求解器混為一談其實 PuLP 只是建模語言真正干活的引擎是 CBC、GLPK、HiGHS 這些求解器。CBC 是 PuLP 默認(rèn)帶的中小規(guī)模問題完全夠用GLPK 更加輕量HiGHS 是目前開源社區(qū)里性能非常強的新銳很多新項目已經(jīng)在用它。如果業(yè)務(wù)規(guī)模到了幾百上千萬變量開源求解器可能力不從心那時才需要考慮 Gurobi 或 CPLEX。但到那個階段之前用 PuLP CBC 已經(jīng)能解決絕大多數(shù)問題。4.4 環(huán)境準(zhǔn)備和最小 Demo安裝很簡單一條命令搞定pip install pulp然后寫下最簡單的線性規(guī)劃程序import pulp prob pulp.LpProblem(demo, pulp.LpMaximize) x pulp.LpVariable(x, lowBound0) y pulp.LpVariable(y, lowBound0) prob 7 * x 5 * y # 目標(biāo)函數(shù) prob 2 * x y 10 # 設(shè)備工時 prob x 3 # A 產(chǎn)品產(chǎn)量上限 prob.solve() print(fx {x.value()}, y {y.value()}) print(fprofit {pulp.value(prob.objective)}) print(pulp.LpStatus[prob.status])運行結(jié)果就是 x3、y4、利潤 41和第 3 節(jié)手算的一致。這個最小 Demo 可以作為以后所有線性規(guī)劃代碼的起點模板。5. 一個完整的排產(chǎn)優(yōu)化案例從建模到落地全流程5.1 業(yè)務(wù)背景與已知條件光說不練沒有用我用一個完整的排產(chǎn)案例把全流程走一遍。某小工廠有兩條產(chǎn)線生產(chǎn) A、B、C 三種產(chǎn)品每種產(chǎn)品的單位利潤、耗用工時、材料用量和需求上限如下表產(chǎn)品產(chǎn)線1工時/件產(chǎn)線2工時/件材料/件利潤/件需求上限A2 小時1 小時5 單位40 元60 件B1 小時2 小時4 單位30 元70 件C1.5 小時1.5 小時3 單位50 元50 件月度可用資源產(chǎn)線1最多 200 小時產(chǎn)線2最多 180 小時原材料最多 500 單位。目標(biāo)是最大化月度總利潤。5.2 建立數(shù)學(xué)模型設(shè) x_A、x_B、x_C 分別代表三種產(chǎn)品的產(chǎn)量模型寫成max Z 40x_A 30x_B 50x_C2x_A x_B 1.5x_C ≤ 200產(chǎn)線1工時x_A 2x_B 1.5x_C ≤ 180產(chǎn)線2工時5x_A 4x_B 3x_C ≤ 500原材料0 ≤ x_A ≤ 600 ≤ x_B ≤ 700 ≤ x_C ≤ 50到這里建模階段就完成了。你會發(fā)現(xiàn)整個過程中我完全沒有去想求解器內(nèi)部怎么迭代只把業(yè)務(wù)規(guī)則翻譯成了公式。5.3 PuLP 完整代碼實現(xiàn)import pulp prob pulp.LpProblem(Monthly_Production_Plan, pulp.LpMaximize) A pulp.LpVariable(A, lowBound0, upBound60, catContinuous) B pulp.LpVariable(B, lowBound0, upBound70, catContinuous) C pulp.LpVariable(C, lowBound0, upBound50, catContinuous) prob 40 * A 30 * B 50 * C, Total_Profit prob 2 * A B 1.5 * C 200, Line1_Capacity prob A 2 * B 1.5 * C 180, Line2_Capacity prob 5 * A 4 * B 3 * C 500, Material_Stock prob.solve() print(status:, pulp.LpStatus[prob.status]) for v in [A, B, C]: print(v.name, , v.value()) print(total profit , pulp.value(prob.objective))運行后輸出status: Optimal A 50.0 B 25.0 C 50.0 total profit 5250.05.4 結(jié)果解讀最優(yōu)解背后的業(yè)務(wù)信息這個結(jié)果非常有意思。C 產(chǎn)品利潤最高直接排滿了上限 50 件A 雖然利潤比 B 高但產(chǎn)量排在 50 件而不是上限 60 件B 只排了 25 件離需求上限 70 件還很遠(yuǎn)。原因在于瓶頸資源。在這個最優(yōu)解里產(chǎn)線1工時和原材料約束都拉滿了分別是 2×50251.5×50200 和 5×504×253×50500而產(chǎn)線2只用了 145 小時剩了 35 小時空閑。這說明真正卡住產(chǎn)能的不是市場需求而是產(chǎn)線1和原材料。進一步還可以看影子價格。我手動算過在這個模型里產(chǎn)線1每增加 1 小時總利潤大約增加 3.33 元原材料每增加 1 單位總利潤大約增加 6.67 元。這就是對偶變量它直接告訴你“哪個瓶頸資源更值得花錢擴充”價值遠(yuǎn)不止算出一個產(chǎn)量方案。6. 求解結(jié)果異常時的排查思路無解、退化與數(shù)值問題6.1 無解的完整排查鏈路求解器報 Infeasible 是新手最常遇到的狀況意思是約束之間互相矛盾可行域是空的。我每次遇到這個報錯都會按下面順序排查檢查所有不等號方向是否寫反尤其是把“至少需要”寫成這種低級錯誤占了無解原因的一大半把約束一條一條注釋掉跑一次看能否恢復(fù)可行這種“二分定位法”通常很快能找到?jīng)_突的那條約束檢查變量上下限是否合理比如某個產(chǎn)線的最大產(chǎn)能寫成了 0那模型當(dāng)然無解檢查量綱是否統(tǒng)一小時和分鐘混用、噸和千克混用都會造成表面上看不出來的沖突。這個排查過程不要急著改代碼先對著模型本身念一遍很多問題在數(shù)學(xué)表達式階段就能發(fā)現(xiàn)。6.2 多重最優(yōu)解目標(biāo)函數(shù)與約束平行有時候模型能解出來但解不止一個。典型場景是目標(biāo)函數(shù)的等值線和某條約束邊界平行這時候所有落在該邊界上的點都是最優(yōu)解利潤完全一樣。對業(yè)務(wù)來說這可能沒問題但也可能很麻煩——比如你有多個同樣利潤的方案但其中某個方案對后續(xù)排班更友好。解決辦法是在目標(biāo)函數(shù)里加一個很小的懲罰項或者偏好項比如在目標(biāo)里減去一個極小量乘以某個變量把求解器“引導(dǎo)”到你更想要的那個解上。這不算作弊這是多目標(biāo)優(yōu)化里很常規(guī)的處理手法。6.3 數(shù)值問題數(shù)量級相差過大的陷阱求解器內(nèi)部用的都是浮點數(shù)對數(shù)值尺度非常敏感。如果模型里同時出現(xiàn) 10? 和 10?? 量級的系數(shù)求解器可能在數(shù)值容差范圍內(nèi)“認(rèn)為”某個約束已經(jīng)滿足但實際誤差很大導(dǎo)致結(jié)果完全失真。處理辦法是統(tǒng)一單位。比如利潤從“元”改成“萬元”工時從“小時”改成“千小時”盡量讓所有系數(shù)落在 0.001 到 1000 這個區(qū)間內(nèi)。還有一種常見情況是懲罰系數(shù)取得太大比如大 M 法里 M 取了 10?結(jié)果把其他約束的精度全沖掉了——M 能取到剛好夠用的程度就行不是越大越好。6.4 規(guī)模變大之后的退路當(dāng)變量和約束數(shù)量漲到幾十萬甚至上百萬單純形法可能會變慢這時候可以考慮 HiGHS 或商業(yè)求解器如果問題還必須加很多整數(shù)變量規(guī)模再大精確求解器也可能撐不住就要考慮列生成、割平面甚至啟發(fā)式算法了。但我的建議是先確保精確模型是對的再談規(guī)模優(yōu)化。很多人一遇到大規(guī)模問題就直接上遺傳算法、模擬退火結(jié)果連最優(yōu)解的參考值都沒有調(diào)參調(diào)到懷疑人生。先用 LP 松弛或者縮小規(guī)模跑一個“理論上的近似最優(yōu)”后面做算法對比才有基準(zhǔn)線。7. 線性規(guī)劃與貪心、動態(tài)規(guī)劃的實際分工7.1 貪心算法特殊結(jié)構(gòu)下的快速通道很多人學(xué)算法時最先接觸貪心像最小生成樹、Dijkstra 最短路這些問題的結(jié)構(gòu)保證了局部最優(yōu)就是全局最優(yōu)。但貪心成立的條件非?？量桃坏┘s束多起來比如同時考慮產(chǎn)能、材料、需求上限你就很難再用“每次都選當(dāng)前最優(yōu)”的思路去湊全局最優(yōu)解。線性規(guī)劃不是要取代貪心而是補上貪心覆蓋不了的那一大片“約束密集”的問題空間。實際問題里貪心適合快速算一個粗略方案線性規(guī)劃適合算精確的最優(yōu)方案。7.2 動態(tài)規(guī)劃狀態(tài)空間能枚舉時的精確解動態(tài)規(guī)劃也很常用核心是把問題拆成重疊子問題用狀態(tài)轉(zhuǎn)移方程逐層求解。只要狀態(tài)空間可控動態(tài)規(guī)劃能給出非常漂亮的精確解。但動態(tài)規(guī)劃的痛點是狀態(tài)爆炸。排產(chǎn)問題如果每個產(chǎn)品的產(chǎn)量范圍有 60 個可能取值、三個產(chǎn)品組合出幾十萬種狀態(tài)再加幾條約束狀態(tài)空間瞬間就收不住了。線性規(guī)劃處理這類問題的思路完全不同變量數(shù)量可以成千上萬求解時間和狀態(tài)枚舉沒有直接關(guān)系。7.3 線性規(guī)劃約束密集、變量連續(xù)時的首選我現(xiàn)在遇到資源分配類需求第一反應(yīng)不是翻算法書找貪心或動態(tài)規(guī)劃而是先問自己三個問題目標(biāo)能不能寫成線性限制能不能寫成線性變量是連續(xù)還是整數(shù)只要前兩個是肯定的線性規(guī)劃就是最穩(wěn)的選擇。除了能給出最優(yōu)解線性規(guī)劃還附送對偶變量、敏感性分析這些額外信息能告訴你“哪個約束是瓶頸”“資源的邊際價值是多少”。這些信息對業(yè)務(wù)決策的價值往往比產(chǎn)量數(shù)字本身更大。7.4 我的選型經(jīng)驗這幾年實踐下來我形成了固定的工作流先花一小時把業(yè)務(wù)翻譯成數(shù)學(xué)模型判斷問題是否是線性能用線性規(guī)劃解決就用線性規(guī)劃變量要求是整數(shù)就升級成 MILP只有問題規(guī)模實在太大、精確求解器扛不住的時候才考慮啟發(fā)式算法。我的切身體會是很多團隊一上來就抱著模擬退火或者遺傳算法不放連最優(yōu)解長什么樣都不知道調(diào)參調(diào)到崩潰。先用線性規(guī)劃跑出一個理論上限再決定要不要上啟發(fā)式這個順序能省下大量無謂的調(diào)參時間。下次再遇到資源分配和排程優(yōu)化的問題我建議你也先試試這個思路。本文還有配套的精品資源點擊獲取