劃與0-1規(guī)劃：Matlab與Lingo建模求解實(shí)戰(zhàn)指南)
1. 從線性到非線性規(guī)劃問題的現(xiàn)實(shí)躍遷在數(shù)學(xué)建模的實(shí)戰(zhàn)中線性規(guī)劃模型因其結(jié)構(gòu)清晰、求解高效往往是我們的首選。無論是經(jīng)典的運(yùn)輸問題、資源分配還是投資組合優(yōu)化線性規(guī)劃都能提供一套漂亮的數(shù)學(xué)框架。然而現(xiàn)實(shí)世界遠(yuǎn)比理想模型復(fù)雜。當(dāng)我們試圖描述成本隨產(chǎn)量非線性增長、利潤與廣告投入呈S型曲線關(guān)系或者決策變量只能取0或1代表“是/否”、“開/關(guān)”時線性規(guī)劃的“直線”世界就立刻顯得捉襟見肘了。這正是“非線性規(guī)劃”與“01規(guī)劃”登場的時刻。它們不是數(shù)學(xué)上的炫技而是解決更真實(shí)、更復(fù)雜決策問題的必要工具。對于參加數(shù)模競賽的同學(xué)而言掌握這兩類模型意味著你的工具箱里多了兩把能處理更棘手問題的“瑞士軍刀”。本文將結(jié)合Matlab和Lingo這兩款利器深入拆解非線性規(guī)劃與01規(guī)劃的核心思想、建模要點(diǎn)與求解策略分享從模型構(gòu)建到代碼實(shí)現(xiàn)的完整經(jīng)驗(yàn)幫你跨越從“知道概念”到“能解出答案”的鴻溝。2. 非線性規(guī)劃當(dāng)目標(biāo)與約束走出“直線”線性規(guī)劃的核心假設(shè)是目標(biāo)函數(shù)和約束條件均為決策變量的線性組合。一旦這個條件被打破我們就進(jìn)入了非線性規(guī)劃的領(lǐng)域。這在實(shí)際問題中極為常見比如經(jīng)濟(jì)學(xué)中的邊際效用遞減目標(biāo)函數(shù)為凹函數(shù)、工程中的最小化阻力目標(biāo)函數(shù)復(fù)雜、化學(xué)反應(yīng)中的平衡濃度約束為非線性方程等。2.1 非線性規(guī)劃模型的標(biāo)準(zhǔn)形式與分類一個標(biāo)準(zhǔn)的非線性規(guī)劃問題可以表述為 求決策變量x使得 最小化或最大化f(x)滿足約束g_i(x) ≤ 0, i 1, ..., m不等式約束h_j(x) 0, j 1, ..., p等式約束x ∈ S通常S為R^n的子集可能包含邊界這里的f(x),g_i(x),h_j(x)至少有一個是非線性函數(shù)。根據(jù)函數(shù)的性質(zhì)非線性規(guī)劃問題可以進(jìn)一步細(xì)分凸規(guī)劃如果f(x)是凸函數(shù)求最小化時且可行域是凸集那么局部最優(yōu)解就是全局最優(yōu)解。這是最“友好”的一類非線性規(guī)劃。非凸規(guī)劃目標(biāo)函數(shù)或可行域非凸。這類問題可能包含多個局部最優(yōu)解找到全局最優(yōu)解非常困難。無約束優(yōu)化只有目標(biāo)函數(shù)沒有約束條件。這是非線性規(guī)劃的基礎(chǔ)。有約束優(yōu)化包含等式或不等式約束求解難度更大。在數(shù)模競賽中我們遇到的大部分是中小規(guī)模、連續(xù)變量的非線性規(guī)劃問題。關(guān)鍵在于如何將實(shí)際問題準(zhǔn)確地轉(zhuǎn)化為這個數(shù)學(xué)形式并選擇合適的工具求解。2.2 Matlab求解非線性規(guī)劃fmincon函數(shù)深度解析Matlab的優(yōu)化工具箱提供了強(qiáng)大的fmincon函數(shù)專門用于求解有約束的非線性多元函數(shù)最小值問題。它的基本調(diào)用格式是[x, fval, exitflag, output] fmincon(fun, x0, A, b, Aeq, beq, lb, ub, nonlcon, options)看起來參數(shù)很多別慌我們逐一拆解其背后的邏輯和實(shí)戰(zhàn)要點(diǎn)fun目標(biāo)函數(shù)句柄。你需要寫一個函數(shù)文件如myfun.m或匿名函數(shù)來計(jì)算f(x)。關(guān)鍵技巧盡量使用向量化運(yùn)算避免在函數(shù)內(nèi)部使用循環(huán)這能極大提升求解速度尤其是在變量維度較高時。x0初始猜測點(diǎn)。這是非線性規(guī)劃求解中最關(guān)鍵也最“玄學(xué)”的一步。因?yàn)閒mincon使用迭代算法通常是內(nèi)點(diǎn)法或序列二次規(guī)劃法不同的初始點(diǎn)可能收斂到不同的局部最優(yōu)解。實(shí)戰(zhàn)經(jīng)驗(yàn)對于非凸問題沒有萬全之策。常用的策略包括1) 根據(jù)問題物理意義給出合理猜測2) 在可行域內(nèi)隨機(jī)生成多個初始點(diǎn)分別求解取最優(yōu)結(jié)果3) 先用全局搜索算法如遺傳算法粗略定位再用fmincon精細(xì)優(yōu)化。A, b, Aeq, beq, lb, ub這些是線性約束和邊界約束。A*x ≤ b和Aeq*x beq定義了線性不等式和等式約束lb和ub是變量的下界和上界。注意即使你的問題包含非線性約束只要存在線性部分也應(yīng)通過這幾個參數(shù)輸入這能幫助求解器更高效地處理問題。nonlcon非線性約束函數(shù)句柄。這個函數(shù)需要返回兩個值不等式約束c(x) ≤ 0和等式約束ceq(x) 0。踩坑提醒務(wù)必確保你的nonlcon函數(shù)能同時計(jì)算c和ceq即使其中一項(xiàng)為空也要返回空數(shù)組[]。函數(shù)定義應(yīng)類似function [c, ceq] mycon(x)。options優(yōu)化選項(xiàng)。這是高手和新手的分水嶺。通過optimoptions(fmincon)可以設(shè)置一系列參數(shù)例如Algorithm選擇算法如interior-point內(nèi)點(diǎn)法默認(rèn)適合大規(guī)模問題、sqp序列二次規(guī)劃適合中小規(guī)模、約束多的問題、active-set有效集法。對于光滑問題interior-point通常是不錯的選擇。Display設(shè)置迭代信息顯示級別iter可以查看每一步的詳細(xì)信息調(diào)試時非常有用。MaxIterations和MaxFunctionEvaluations防止程序陷入無限循環(huán)或計(jì)算時間過長。OptimalityTolerance和StepTolerance收斂容差。如果結(jié)果精度不夠可以適當(dāng)調(diào)小這些值如1e-8。一個完整的建模與求解示例假設(shè)我們要優(yōu)化一個產(chǎn)品生產(chǎn)問題利潤函數(shù)為f(x) - (2*x1 3*x2 x1*x2)求最大利潤即求-f的最小值受限于資源約束x1^2 x2^2 ≤ 4和非負(fù)條件x1, x2 ≥ 0。% 1. 定義目標(biāo)函數(shù) (求最小化所以是 -利潤) fun (x) - (2*x(1) 3*x(2) x(1)*x(2)); % 2. 定義非線性約束 x1^2 x2^2 - 4 0 function [c, ceq] circlecon(x) c x(1)^2 x(2)^2 - 4; % 不等式約束要求 c 0 ceq []; % 沒有等式約束 end % 3. 設(shè)置其他參數(shù) x0 [1, 1]; % 初始猜測 A []; b []; Aeq []; beq []; % 無線性約束 lb [0, 0]; % 下界 ub []; % 無上界 % 4. 調(diào)用 fmincon 求解 options optimoptions(fmincon, Display, iter, Algorithm, interior-point); [x_opt, fval_opt] fmincon(fun, x0, A, b, Aeq, beq, lb, ub, circlecon, options); % 5. 輸出結(jié)果 fprintf(最優(yōu)解: x1 %.4f, x2 %.4f\n, x_opt(1), x_opt(2)); fprintf(最大利潤為: %.4f\n, -fval_opt); % 注意取負(fù)號運(yùn)行這段代碼你會看到迭代過程并得到最優(yōu)解。通過調(diào)整x0為[0,0]或[2,0]你可以觀察是否收斂到同一個點(diǎn)以此初步判斷問題的凸性。2.3 Lingo求解非線性規(guī)劃更貼近自然語言的建模Lingo的魅力在于其建模語言幾乎是對數(shù)學(xué)模型的直接翻譯特別適合快速原型驗(yàn)證。對于上面的例子在Lingo中的模型文件.lg4可以這樣寫MODEL: ! 定義集合和變量; SETS: PRODUCT /1..2/: X, ProfitCoeff; ENDSETS DATA: ProfitCoeff 2, 3; ! 利潤系數(shù); ENDDATA ! 目標(biāo)函數(shù)最大化總利潤; MAX SUM(PRODUCT(i): ProfitCoeff(i) * X(i)) X(1)*X(2); ! 資源約束; X(1)^2 X(2)^2 4; ! 非負(fù)約束; FOR(PRODUCT(i): X(i) 0); END點(diǎn)擊求解Lingo會調(diào)用其非線性求解器通常是廣義既約梯度法進(jìn)行計(jì)算。Lingo的優(yōu)勢在于1) 語法直觀易于檢查和修改模型2) 對于中小規(guī)模問題設(shè)置簡單3) 結(jié)果報告清晰包含靈敏度分析對于線性部分。但需要注意Lingo的非線性全局求解能力有限對于非凸問題也可能只找到局部最優(yōu)解。在Lingo中可以通過LINGO - Options - Global Solver勾選全局求解器來嘗試尋找全局最優(yōu)但這會顯著增加計(jì)算時間。Matlab與Lingo的選擇心得如果你的問題需要集成到更大的算法流程中如與仿真、數(shù)據(jù)處理結(jié)合或者需要高度定制化的求解流程和算法Matlab是不二之選。如果你的核心工作是快速、清晰地建立并求解一個獨(dú)立的優(yōu)化模型特別是向非編程背景的隊(duì)友或評委展示模型時Lingo的代碼可讀性更具優(yōu)勢。在數(shù)模比賽中可以根據(jù)團(tuán)隊(duì)技能和問題特點(diǎn)靈活選用甚至用Lingo快速驗(yàn)證模型正確性再用Matlab進(jìn)行深入分析或集成。3. 01規(guī)劃離散決策的利器當(dāng)決策變量不是連續(xù)的而是只能取0或1時我們就進(jìn)入了整數(shù)規(guī)劃的特例——01規(guī)劃Binary Programming的領(lǐng)域。這用來表示一系列“是或否”、“選擇或不選擇”、“打開或關(guān)閉”的決策。例如選址問題某個地點(diǎn)建或不建工廠、投資選擇某個項(xiàng)目投或不投、背包問題某件物品帶或不帶、人員排班某個時段是否安排某人上班等。3.1 01規(guī)劃模型的特點(diǎn)與挑戰(zhàn)01規(guī)劃模型在形式上與線性或非線性規(guī)劃類似只是增加了x_i ∈ {0, 1}的約束。正是這個簡單的約束將問題從“多項(xiàng)式時間可解”的領(lǐng)域?qū)τ诰€性規(guī)劃拖入了“NP-Hard”的復(fù)雜世界。求解01規(guī)劃的核心挑戰(zhàn)在于組合爆炸n個01變量會產(chǎn)生2^n種可能的組合。窮舉法對于稍大的n就完全不可行。因此求解01規(guī)劃主要依賴兩類方法精確算法如分支定界法、割平面法。這類方法能保證找到全局最優(yōu)解但最壞情況下的計(jì)算時間仍可能很長。啟發(fā)式/元啟發(fā)式算法如遺傳算法、模擬退火、禁忌搜索。這類方法不能保證找到全局最優(yōu)但能在可接受的時間內(nèi)找到高質(zhì)量接近最優(yōu)的解非常適合大規(guī)模或復(fù)雜的01規(guī)劃問題。在數(shù)模競賽中我們通常借助工具內(nèi)置的求解器它們封裝了成熟的精確算法。3.2 Matlab求解01規(guī)劃intlinprog與ga的配合對于線性的01規(guī)劃Matlab的優(yōu)化工具箱提供了intlinprog函數(shù)。它是linprog的整數(shù)規(guī)劃版本可以高效求解混合整數(shù)線性規(guī)劃問題其中就包含所有變量都是0-1的情況。基本調(diào)用格式[x, fval] intlinprog(f, intcon, A, b, Aeq, beq, lb, ub)f目標(biāo)函數(shù)系數(shù)向量。intcon指定哪些變量是整數(shù)變量。對于純01規(guī)劃intcon就是所有變量的索引例如1:n。A, b, Aeq, beq, lb, ub線性約束和邊界。這里有個關(guān)鍵點(diǎn)為了將變量限制在0和1我們必須同時設(shè)置lb zeros(n,1)和ub ones(n,1)。intlinprog會結(jié)合這些邊界和intcon參數(shù)將變量識別為01變量。示例經(jīng)典的背包問題。有5件物品重量w[2,3,4,5,9]價值v[3,4,5,8,10]背包容量C15。選擇哪些物品使得總價值最大且總重量不超過容量f -v; % 求最大價值轉(zhuǎn)化為求最小化 -v*x intcon 1:5; A w; b C; Aeq []; beq []; lb zeros(5,1); ub ones(5,1); [x_opt, fval_opt] intlinprog(f, intcon, A, b, Aeq, beq, lb, ub); disp(選擇的物品索引); find(x_opt 0.5) % 由于數(shù)值計(jì)算解可能接近0或1但不完全等于 disp([最大價值, num2str(-fval_opt)]);對于非線性的01規(guī)劃即目標(biāo)函數(shù)或約束包含非線性項(xiàng)且變量為01情況就復(fù)雜多了。Matlab沒有專門的函數(shù)。一個常見的處理思路是使用遺傳算法。雖然遺傳算法不要求問題可微或連續(xù)但它處理嚴(yán)格的01約束和復(fù)雜非線性約束的能力更強(qiáng)。使用ga求解非線性01規(guī)劃示例假設(shè)一個簡單的非線性01規(guī)劃max x1*x2 x3滿足x1 2*x2 - x3 1x1, x2, x3 ∈ {0,1}。% 定義適應(yīng)度函數(shù)求最大所以取負(fù) fun (x) - (x(1)*x(2) x(3)); % 變量個數(shù) nvars 3; % 線性不等式約束 A*x b A [1, 2, -1]; b 1; % 定義變量的上下界為0和1 lb [0, 0, 0]; ub [1, 1, 1]; % 關(guān)鍵使用自定義的整數(shù)約束創(chuàng)建函數(shù) function [state, options, optchanged] binaryconstraint(options, state, flag) optchanged false; if strcmp(flag, iter) % 將種群中所有個體的變量舍入到最接近的0或1 state.Population round(state.Population); end end % 設(shè)置遺傳算法選項(xiàng)加入自定義輸出函數(shù)來強(qiáng)制01約束 options optimoptions(ga, ... Display, iter, ... ConstraintTolerance, 1e-6, ... PlotFcn, gaplotbestf, ... OutputFcn, binaryconstraint); % 關(guān)鍵加入輸出函數(shù) % 調(diào)用ga求解。注意ga默認(rèn)處理邊界約束但線性約束A,b也需要傳入。 [x_opt, fval_opt] ga(fun, nvars, A, b, [], [], lb, ub, [], [], options); fprintf(最優(yōu)解: [%d, %d, %d]\n, round(x_opt)); fprintf(最優(yōu)值: %.4f\n, -fval_opt);重要提示這種方法在迭代中強(qiáng)制舍入是一種啟發(fā)式處理它破壞了遺傳算法的自然進(jìn)化過程可能影響找到全局最優(yōu)解的能力并且不能嚴(yán)格保證線性約束在舍入后仍然滿足。對于復(fù)雜的非線性01規(guī)劃這通常是一個折衷方案。更嚴(yán)謹(jǐn)?shù)淖龇ㄊ窃O(shè)計(jì)特殊的編碼方式和遺傳算子來保證01屬性。3.3 Lingo求解01規(guī)劃語法簡潔直擊核心在Lingo中處理01規(guī)劃非常直接只需在變量定義后加上BIN函數(shù)即可。以上述非線性01規(guī)劃為例MODEL: SETS: ITEM /1..3/: X; ENDSETS ! 目標(biāo)函數(shù); MAX X(1) * X(2) X(3); ! 線性約束; X(1) 2*X(2) - X(3) 1; ! 01變量聲明; FOR(ITEM(i): BIN(X(i))); END點(diǎn)擊求解Lingo會調(diào)用其整數(shù)規(guī)劃求解器通?；诜种Фń绶ㄟM(jìn)行求解。對于非線性01規(guī)劃Lingo會先嘗試線性化如果可能或者使用其全局求解器。Lingo的優(yōu)勢再次凸顯建模極其簡潔BIN一句聲明就搞定省去了在Matlab中處理邊界和整數(shù)約束的麻煩。對于混合整數(shù)非線性規(guī)劃Lingo的求解能力往往比Matlab的內(nèi)置函數(shù)更穩(wěn)健和方便。01規(guī)劃建模的實(shí)用技巧邏輯約束的轉(zhuǎn)化很多邏輯關(guān)系可以用01變量和線性約束來表達(dá)。例如“如果項(xiàng)目A被選中(x_A1)則項(xiàng)目B也必須被選中(x_B1)”x_A x_B。“在項(xiàng)目A和B中至少選一個”x_A x_B 1?！霸陧?xiàng)目A和B中至多選一個”x_A x_B 1?！绊?xiàng)目C當(dāng)且僅當(dāng)項(xiàng)目A和B都選中時才被選中”2*x_C x_A x_B且x_A x_B - 1 x_C。 熟練掌握這些轉(zhuǎn)化是建立復(fù)雜01規(guī)劃模型的基本功。Big-M法用于處理帶有固定成本的決策或者將非線性關(guān)系如if-then線性化。例如如果選擇生產(chǎn)某種產(chǎn)品y1會產(chǎn)生固定成本F且產(chǎn)量x有上限M。則可以寫成x M * y并且目標(biāo)函數(shù)中包含F(xiàn) * y。這里M是一個足夠大的數(shù)當(dāng)y0時強(qiáng)制x0當(dāng)y1時x可以取到上限M以內(nèi)的任何值。4. 數(shù)模實(shí)戰(zhàn)融合非線性與01規(guī)劃的綜合應(yīng)用與排錯在真正的數(shù)模賽題中純非線性或純01規(guī)劃的問題較少更多是兩者的混合或者與其他模型如動態(tài)規(guī)劃、圖論結(jié)合。例如一個設(shè)施選址問題01決策中每個設(shè)施的運(yùn)營成本可能是其服務(wù)量的非線性函數(shù)。4.1 典型賽題思路拆解假設(shè)一個簡化版的“電動汽車充電站選址與容量規(guī)劃”問題01決策在若干個候選位置中選擇哪些位置建設(shè)充電站y_j ∈ {0,1}。連續(xù)決策每個充電站的容量充電樁數(shù)量x_j連續(xù)變量。非線性關(guān)系建設(shè)成本可能與容量呈規(guī)模經(jīng)濟(jì)效應(yīng)凹函數(shù)如cost_j a * sqrt(x_j) b或者擁堵成本凸函數(shù)。約束滿足所有區(qū)域的需求容量總和有上限投資總預(yù)算限制等。建模步驟定義集合候選站址集合J需求區(qū)域集合I。定義變量01變量y_j連續(xù)變量x_j以及可能的需求分配變量z_ij從站j滿足區(qū)域i的需求量。目標(biāo)函數(shù)最小化總成本 總建設(shè)成本非線性含y_j和x_j 總運(yùn)營/輸電成本可能是z_ij的函數(shù)。約束需求滿足對每個區(qū)域i∑_j z_ij demand_i。容量限制對每個站j∑_i z_ij x_j且x_j M * y_jBig-M法如果y_j0則x_j0。邏輯約束例如某個區(qū)域必須被至少一個站覆蓋∑_j a_ij * y_j 1其中a_ij表示站j是否能覆蓋區(qū)域i。預(yù)算約束總建設(shè)成本∑_j (F_j * y_j f(x_j)) ≤ Budget。求解策略這是一個混合整數(shù)非線性規(guī)劃問題。如果非線性部分可以線性化或分段線性化可以嘗試用Lingo或Matlab的intlinprog結(jié)合線性化技巧。如果非線性部分復(fù)雜可以考慮用啟發(fā)式算法如遺傳算法同時優(yōu)化y和x或者在y固定的情況下x的子問題是一個連續(xù)非線性規(guī)劃可以交替優(yōu)化。4.2 常見錯誤與調(diào)試心得在實(shí)現(xiàn)和求解這類模型時新手常會踩一些坑“無可行解”錯誤這是最令人頭疼的。首先檢查所有約束是否自相矛盾。例如邊界lb ub或者兩個約束聯(lián)合起來使得可行域?yàn)榭?。調(diào)試方法逐步注釋掉部分約束特別是非線性約束和復(fù)雜的邏輯約束先讓模型有解再逐個加入約束定位問題源。在Lingo中可以使用LINGO - Generate - Display model查看完整的線性化后的模型檢查約束。在Matlab中檢查A,b,Aeq,beq,lb,ub的維度是否正確?！敖獾馁|(zhì)量差”或“陷入局部最優(yōu)”對于非線性規(guī)劃這通常與初始點(diǎn)x0有關(guān)。策略進(jìn)行多初始點(diǎn)搜索。寫一個循環(huán)隨機(jī)生成多個x0確保在邊界內(nèi)或滿足簡單約束分別調(diào)用fmincon記錄最優(yōu)解。對于01規(guī)劃如果使用啟發(fā)式算法可以增加種群大小和迭代次數(shù)。Lingo報錯“NLP Solver failed”或求解時間過長非線性問題可能非凸Lingo的默認(rèn)本地求解器卡住了。嘗試在Lingo菜單LINGO - Options - Global Solver中勾選Use Global Solver。這會啟用全局優(yōu)化但計(jì)算時間會大幅增加。對于大規(guī)模問題這可能不現(xiàn)實(shí)此時需要考慮問題重構(gòu)或采用啟發(fā)式方法。數(shù)值不穩(wěn)定目標(biāo)函數(shù)或約束條件尺度差異巨大例如一個變量范圍是[0, 1]另一個是[10000, 20000]會導(dǎo)致求解器數(shù)值計(jì)算困難收斂緩慢或不準(zhǔn)確。解決方案對變量進(jìn)行縮放使其處于相近的數(shù)量級上例如都縮放到[0, 1]或[-1, 1]附近。這在Matlab和Lingo中都同樣重要。模型正確但求解器不收斂檢查優(yōu)化選項(xiàng)。在Matlab中適當(dāng)增加MaxIterations和MaxFunctionEvaluations。檢查收斂容差OptimalityTolerance和StepTolerance如果設(shè)置得太小可能永遠(yuǎn)達(dá)不到。有時候稍微放松容差如從1e-10調(diào)到1e-6就能讓求解器成功終止并獲得一個可接受的解。最后分享一個個人在數(shù)模競賽中處理復(fù)雜規(guī)劃問題的習(xí)慣永遠(yuǎn)從最簡單、最核心的模型版本開始。先忽略次要的非線性項(xiàng)用線性模型和01變量把主體邏輯跑通得到基準(zhǔn)解和運(yùn)行時間。然后再逐步加入非線性部分、更復(fù)雜的約束并觀察解的變化和計(jì)算時間的增長。這樣既能保證在有限時間內(nèi)有一個保底的模型也能有條理地評估模型復(fù)雜化帶來的收益與成本。記住在數(shù)模比賽中一個能跑出合理結(jié)果、邏輯清晰的簡化模型遠(yuǎn)勝過一個理論上完美但無法求解或求解不穩(wěn)定的復(fù)雜模型。 非線性規(guī)劃與0-1規(guī)劃：Matlab與Lingo建模求解實(shí)戰(zhàn)指南
1. 從線性到非線性規(guī)劃問題的現(xiàn)實(shí)躍遷在數(shù)學(xué)建模的實(shí)戰(zhàn)中線性規(guī)劃模型因其結(jié)構(gòu)清晰、求解高效往往是我們的首選。無論是經(jīng)典的運(yùn)輸問題、資源分配還是投資組合優(yōu)化線性規(guī)劃都能提供一套漂亮的數(shù)學(xué)框架。然而現(xiàn)實(shí)世界遠(yuǎn)比理想模型復(fù)雜。當(dāng)我們試圖描述成本隨產(chǎn)量非線性增長、利潤與廣告投入呈S型曲線關(guān)系或者決策變量只能取0或1代表“是/否”、“開/關(guān)”時線性規(guī)劃的“直線”世界就立刻顯得捉襟見肘了。這正是“非線性規(guī)劃”與“01規(guī)劃”登場的時刻。它們不是數(shù)學(xué)上的炫技而是解決更真實(shí)、更復(fù)雜決策問題的必要工具。對于參加數(shù)模競賽的同學(xué)而言掌握這兩類模型意味著你的工具箱里多了兩把能處理更棘手問題的“瑞士軍刀”。本文將結(jié)合Matlab和Lingo這兩款利器深入拆解非線性規(guī)劃與01規(guī)劃的核心思想、建模要點(diǎn)與求解策略分享從模型構(gòu)建到代碼實(shí)現(xiàn)的完整經(jīng)驗(yàn)幫你跨越從“知道概念”到“能解出答案”的鴻溝。2. 非線性規(guī)劃當(dāng)目標(biāo)與約束走出“直線”線性規(guī)劃的核心假設(shè)是目標(biāo)函數(shù)和約束條件均為決策變量的線性組合。一旦這個條件被打破我們就進(jìn)入了非線性規(guī)劃的領(lǐng)域。這在實(shí)際問題中極為常見比如經(jīng)濟(jì)學(xué)中的邊際效用遞減目標(biāo)函數(shù)為凹函數(shù)、工程中的最小化阻力目標(biāo)函數(shù)復(fù)雜、化學(xué)反應(yīng)中的平衡濃度約束為非線性方程等。2.1 非線性規(guī)劃模型的標(biāo)準(zhǔn)形式與分類一個標(biāo)準(zhǔn)的非線性規(guī)劃問題可以表述為 求決策變量x使得 最小化或最大化f(x)滿足約束g_i(x) ≤ 0, i 1, ..., m不等式約束h_j(x) 0, j 1, ..., p等式約束x ∈ S通常S為R^n的子集可能包含邊界這里的f(x),g_i(x),h_j(x)至少有一個是非線性函數(shù)。根據(jù)函數(shù)的性質(zhì)非線性規(guī)劃問題可以進(jìn)一步細(xì)分凸規(guī)劃如果f(x)是凸函數(shù)求最小化時且可行域是凸集那么局部最優(yōu)解就是全局最優(yōu)解。這是最“友好”的一類非線性規(guī)劃。非凸規(guī)劃目標(biāo)函數(shù)或可行域非凸。這類問題可能包含多個局部最優(yōu)解找到全局最優(yōu)解非常困難。無約束優(yōu)化只有目標(biāo)函數(shù)沒有約束條件。這是非線性規(guī)劃的基礎(chǔ)。有約束優(yōu)化包含等式或不等式約束求解難度更大。在數(shù)模競賽中我們遇到的大部分是中小規(guī)模、連續(xù)變量的非線性規(guī)劃問題。關(guān)鍵在于如何將實(shí)際問題準(zhǔn)確地轉(zhuǎn)化為這個數(shù)學(xué)形式并選擇合適的工具求解。2.2 Matlab求解非線性規(guī)劃fmincon函數(shù)深度解析Matlab的優(yōu)化工具箱提供了強(qiáng)大的fmincon函數(shù)專門用于求解有約束的非線性多元函數(shù)最小值問題。它的基本調(diào)用格式是[x, fval, exitflag, output] fmincon(fun, x0, A, b, Aeq, beq, lb, ub, nonlcon, options)看起來參數(shù)很多別慌我們逐一拆解其背后的邏輯和實(shí)戰(zhàn)要點(diǎn)fun目標(biāo)函數(shù)句柄。你需要寫一個函數(shù)文件如myfun.m或匿名函數(shù)來計(jì)算f(x)。關(guān)鍵技巧盡量使用向量化運(yùn)算避免在函數(shù)內(nèi)部使用循環(huán)這能極大提升求解速度尤其是在變量維度較高時。x0初始猜測點(diǎn)。這是非線性規(guī)劃求解中最關(guān)鍵也最“玄學(xué)”的一步。因?yàn)閒mincon使用迭代算法通常是內(nèi)點(diǎn)法或序列二次規(guī)劃法不同的初始點(diǎn)可能收斂到不同的局部最優(yōu)解。實(shí)戰(zhàn)經(jīng)驗(yàn)對于非凸問題沒有萬全之策。常用的策略包括1) 根據(jù)問題物理意義給出合理猜測2) 在可行域內(nèi)隨機(jī)生成多個初始點(diǎn)分別求解取最優(yōu)結(jié)果3) 先用全局搜索算法如遺傳算法粗略定位再用fmincon精細(xì)優(yōu)化。A, b, Aeq, beq, lb, ub這些是線性約束和邊界約束。A*x ≤ b和Aeq*x beq定義了線性不等式和等式約束lb和ub是變量的下界和上界。注意即使你的問題包含非線性約束只要存在線性部分也應(yīng)通過這幾個參數(shù)輸入這能幫助求解器更高效地處理問題。nonlcon非線性約束函數(shù)句柄。這個函數(shù)需要返回兩個值不等式約束c(x) ≤ 0和等式約束ceq(x) 0。踩坑提醒務(wù)必確保你的nonlcon函數(shù)能同時計(jì)算c和ceq即使其中一項(xiàng)為空也要返回空數(shù)組[]。函數(shù)定義應(yīng)類似function [c, ceq] mycon(x)。options優(yōu)化選項(xiàng)。這是高手和新手的分水嶺。通過optimoptions(fmincon)可以設(shè)置一系列參數(shù)例如Algorithm選擇算法如interior-point內(nèi)點(diǎn)法默認(rèn)適合大規(guī)模問題、sqp序列二次規(guī)劃適合中小規(guī)模、約束多的問題、active-set有效集法。對于光滑問題interior-point通常是不錯的選擇。Display設(shè)置迭代信息顯示級別iter可以查看每一步的詳細(xì)信息調(diào)試時非常有用。MaxIterations和MaxFunctionEvaluations防止程序陷入無限循環(huán)或計(jì)算時間過長。OptimalityTolerance和StepTolerance收斂容差。如果結(jié)果精度不夠可以適當(dāng)調(diào)小這些值如1e-8。一個完整的建模與求解示例假設(shè)我們要優(yōu)化一個產(chǎn)品生產(chǎn)問題利潤函數(shù)為f(x) - (2*x1 3*x2 x1*x2)求最大利潤即求-f的最小值受限于資源約束x1^2 x2^2 ≤ 4和非負(fù)條件x1, x2 ≥ 0。% 1. 定義目標(biāo)函數(shù) (求最小化所以是 -利潤) fun (x) - (2*x(1) 3*x(2) x(1)*x(2)); % 2. 定義非線性約束 x1^2 x2^2 - 4 0 function [c, ceq] circlecon(x) c x(1)^2 x(2)^2 - 4; % 不等式約束要求 c 0 ceq []; % 沒有等式約束 end % 3. 設(shè)置其他參數(shù) x0 [1, 1]; % 初始猜測 A []; b []; Aeq []; beq []; % 無線性約束 lb [0, 0]; % 下界 ub []; % 無上界 % 4. 調(diào)用 fmincon 求解 options optimoptions(fmincon, Display, iter, Algorithm, interior-point); [x_opt, fval_opt] fmincon(fun, x0, A, b, Aeq, beq, lb, ub, circlecon, options); % 5. 輸出結(jié)果 fprintf(最優(yōu)解: x1 %.4f, x2 %.4f\n, x_opt(1), x_opt(2)); fprintf(最大利潤為: %.4f\n, -fval_opt); % 注意取負(fù)號運(yùn)行這段代碼你會看到迭代過程并得到最優(yōu)解。通過調(diào)整x0為[0,0]或[2,0]你可以觀察是否收斂到同一個點(diǎn)以此初步判斷問題的凸性。2.3 Lingo求解非線性規(guī)劃更貼近自然語言的建模Lingo的魅力在于其建模語言幾乎是對數(shù)學(xué)模型的直接翻譯特別適合快速原型驗(yàn)證。對于上面的例子在Lingo中的模型文件.lg4可以這樣寫MODEL: ! 定義集合和變量; SETS: PRODUCT /1..2/: X, ProfitCoeff; ENDSETS DATA: ProfitCoeff 2, 3; ! 利潤系數(shù); ENDDATA ! 目標(biāo)函數(shù)最大化總利潤; MAX SUM(PRODUCT(i): ProfitCoeff(i) * X(i)) X(1)*X(2); ! 資源約束; X(1)^2 X(2)^2 4; ! 非負(fù)約束; FOR(PRODUCT(i): X(i) 0); END點(diǎn)擊求解Lingo會調(diào)用其非線性求解器通常是廣義既約梯度法進(jìn)行計(jì)算。Lingo的優(yōu)勢在于1) 語法直觀易于檢查和修改模型2) 對于中小規(guī)模問題設(shè)置簡單3) 結(jié)果報告清晰包含靈敏度分析對于線性部分。但需要注意Lingo的非線性全局求解能力有限對于非凸問題也可能只找到局部最優(yōu)解。在Lingo中可以通過LINGO - Options - Global Solver勾選全局求解器來嘗試尋找全局最優(yōu)但這會顯著增加計(jì)算時間。Matlab與Lingo的選擇心得如果你的問題需要集成到更大的算法流程中如與仿真、數(shù)據(jù)處理結(jié)合或者需要高度定制化的求解流程和算法Matlab是不二之選。如果你的核心工作是快速、清晰地建立并求解一個獨(dú)立的優(yōu)化模型特別是向非編程背景的隊(duì)友或評委展示模型時Lingo的代碼可讀性更具優(yōu)勢。在數(shù)模比賽中可以根據(jù)團(tuán)隊(duì)技能和問題特點(diǎn)靈活選用甚至用Lingo快速驗(yàn)證模型正確性再用Matlab進(jìn)行深入分析或集成。3. 01規(guī)劃離散決策的利器當(dāng)決策變量不是連續(xù)的而是只能取0或1時我們就進(jìn)入了整數(shù)規(guī)劃的特例——01規(guī)劃Binary Programming的領(lǐng)域。這用來表示一系列“是或否”、“選擇或不選擇”、“打開或關(guān)閉”的決策。例如選址問題某個地點(diǎn)建或不建工廠、投資選擇某個項(xiàng)目投或不投、背包問題某件物品帶或不帶、人員排班某個時段是否安排某人上班等。3.1 01規(guī)劃模型的特點(diǎn)與挑戰(zhàn)01規(guī)劃模型在形式上與線性或非線性規(guī)劃類似只是增加了x_i ∈ {0, 1}的約束。正是這個簡單的約束將問題從“多項(xiàng)式時間可解”的領(lǐng)域?qū)τ诰€性規(guī)劃拖入了“NP-Hard”的復(fù)雜世界。求解01規(guī)劃的核心挑戰(zhàn)在于組合爆炸n個01變量會產(chǎn)生2^n種可能的組合。窮舉法對于稍大的n就完全不可行。因此求解01規(guī)劃主要依賴兩類方法精確算法如分支定界法、割平面法。這類方法能保證找到全局最優(yōu)解但最壞情況下的計(jì)算時間仍可能很長。啟發(fā)式/元啟發(fā)式算法如遺傳算法、模擬退火、禁忌搜索。這類方法不能保證找到全局最優(yōu)但能在可接受的時間內(nèi)找到高質(zhì)量接近最優(yōu)的解非常適合大規(guī)模或復(fù)雜的01規(guī)劃問題。在數(shù)模競賽中我們通常借助工具內(nèi)置的求解器它們封裝了成熟的精確算法。3.2 Matlab求解01規(guī)劃intlinprog與ga的配合對于線性的01規(guī)劃Matlab的優(yōu)化工具箱提供了intlinprog函數(shù)。它是linprog的整數(shù)規(guī)劃版本可以高效求解混合整數(shù)線性規(guī)劃問題其中就包含所有變量都是0-1的情況。基本調(diào)用格式[x, fval] intlinprog(f, intcon, A, b, Aeq, beq, lb, ub)f目標(biāo)函數(shù)系數(shù)向量。intcon指定哪些變量是整數(shù)變量。對于純01規(guī)劃intcon就是所有變量的索引例如1:n。A, b, Aeq, beq, lb, ub線性約束和邊界。這里有個關(guān)鍵點(diǎn)為了將變量限制在0和1我們必須同時設(shè)置lb zeros(n,1)和ub ones(n,1)。intlinprog會結(jié)合這些邊界和intcon參數(shù)將變量識別為01變量。示例經(jīng)典的背包問題。有5件物品重量w[2,3,4,5,9]價值v[3,4,5,8,10]背包容量C15。選擇哪些物品使得總價值最大且總重量不超過容量f -v; % 求最大價值轉(zhuǎn)化為求最小化 -v*x intcon 1:5; A w; b C; Aeq []; beq []; lb zeros(5,1); ub ones(5,1); [x_opt, fval_opt] intlinprog(f, intcon, A, b, Aeq, beq, lb, ub); disp(選擇的物品索引); find(x_opt 0.5) % 由于數(shù)值計(jì)算解可能接近0或1但不完全等于 disp([最大價值, num2str(-fval_opt)]);對于非線性的01規(guī)劃即目標(biāo)函數(shù)或約束包含非線性項(xiàng)且變量為01情況就復(fù)雜多了。Matlab沒有專門的函數(shù)。一個常見的處理思路是使用遺傳算法。雖然遺傳算法不要求問題可微或連續(xù)但它處理嚴(yán)格的01約束和復(fù)雜非線性約束的能力更強(qiáng)。使用ga求解非線性01規(guī)劃示例假設(shè)一個簡單的非線性01規(guī)劃max x1*x2 x3滿足x1 2*x2 - x3 1x1, x2, x3 ∈ {0,1}。% 定義適應(yīng)度函數(shù)求最大所以取負(fù) fun (x) - (x(1)*x(2) x(3)); % 變量個數(shù) nvars 3; % 線性不等式約束 A*x b A [1, 2, -1]; b 1; % 定義變量的上下界為0和1 lb [0, 0, 0]; ub [1, 1, 1]; % 關(guān)鍵使用自定義的整數(shù)約束創(chuàng)建函數(shù) function [state, options, optchanged] binaryconstraint(options, state, flag) optchanged false; if strcmp(flag, iter) % 將種群中所有個體的變量舍入到最接近的0或1 state.Population round(state.Population); end end % 設(shè)置遺傳算法選項(xiàng)加入自定義輸出函數(shù)來強(qiáng)制01約束 options optimoptions(ga, ... Display, iter, ... ConstraintTolerance, 1e-6, ... PlotFcn, gaplotbestf, ... OutputFcn, binaryconstraint); % 關(guān)鍵加入輸出函數(shù) % 調(diào)用ga求解。注意ga默認(rèn)處理邊界約束但線性約束A,b也需要傳入。 [x_opt, fval_opt] ga(fun, nvars, A, b, [], [], lb, ub, [], [], options); fprintf(最優(yōu)解: [%d, %d, %d]\n, round(x_opt)); fprintf(最優(yōu)值: %.4f\n, -fval_opt);重要提示這種方法在迭代中強(qiáng)制舍入是一種啟發(fā)式處理它破壞了遺傳算法的自然進(jìn)化過程可能影響找到全局最優(yōu)解的能力并且不能嚴(yán)格保證線性約束在舍入后仍然滿足。對于復(fù)雜的非線性01規(guī)劃這通常是一個折衷方案。更嚴(yán)謹(jǐn)?shù)淖龇ㄊ窃O(shè)計(jì)特殊的編碼方式和遺傳算子來保證01屬性。3.3 Lingo求解01規(guī)劃語法簡潔直擊核心在Lingo中處理01規(guī)劃非常直接只需在變量定義后加上BIN函數(shù)即可。以上述非線性01規(guī)劃為例MODEL: SETS: ITEM /1..3/: X; ENDSETS ! 目標(biāo)函數(shù); MAX X(1) * X(2) X(3); ! 線性約束; X(1) 2*X(2) - X(3) 1; ! 01變量聲明; FOR(ITEM(i): BIN(X(i))); END點(diǎn)擊求解Lingo會調(diào)用其整數(shù)規(guī)劃求解器通?；诜种Фń绶ㄟM(jìn)行求解。對于非線性01規(guī)劃Lingo會先嘗試線性化如果可能或者使用其全局求解器。Lingo的優(yōu)勢再次凸顯建模極其簡潔BIN一句聲明就搞定省去了在Matlab中處理邊界和整數(shù)約束的麻煩。對于混合整數(shù)非線性規(guī)劃Lingo的求解能力往往比Matlab的內(nèi)置函數(shù)更穩(wěn)健和方便。01規(guī)劃建模的實(shí)用技巧邏輯約束的轉(zhuǎn)化很多邏輯關(guān)系可以用01變量和線性約束來表達(dá)。例如“如果項(xiàng)目A被選中(x_A1)則項(xiàng)目B也必須被選中(x_B1)”x_A x_B。“在項(xiàng)目A和B中至少選一個”x_A x_B 1?！霸陧?xiàng)目A和B中至多選一個”x_A x_B 1?！绊?xiàng)目C當(dāng)且僅當(dāng)項(xiàng)目A和B都選中時才被選中”2*x_C x_A x_B且x_A x_B - 1 x_C。 熟練掌握這些轉(zhuǎn)化是建立復(fù)雜01規(guī)劃模型的基本功。Big-M法用于處理帶有固定成本的決策或者將非線性關(guān)系如if-then線性化。例如如果選擇生產(chǎn)某種產(chǎn)品y1會產(chǎn)生固定成本F且產(chǎn)量x有上限M。則可以寫成x M * y并且目標(biāo)函數(shù)中包含F(xiàn) * y。這里M是一個足夠大的數(shù)當(dāng)y0時強(qiáng)制x0當(dāng)y1時x可以取到上限M以內(nèi)的任何值。4. 數(shù)模實(shí)戰(zhàn)融合非線性與01規(guī)劃的綜合應(yīng)用與排錯在真正的數(shù)模賽題中純非線性或純01規(guī)劃的問題較少更多是兩者的混合或者與其他模型如動態(tài)規(guī)劃、圖論結(jié)合。例如一個設(shè)施選址問題01決策中每個設(shè)施的運(yùn)營成本可能是其服務(wù)量的非線性函數(shù)。4.1 典型賽題思路拆解假設(shè)一個簡化版的“電動汽車充電站選址與容量規(guī)劃”問題01決策在若干個候選位置中選擇哪些位置建設(shè)充電站y_j ∈ {0,1}。連續(xù)決策每個充電站的容量充電樁數(shù)量x_j連續(xù)變量。非線性關(guān)系建設(shè)成本可能與容量呈規(guī)模經(jīng)濟(jì)效應(yīng)凹函數(shù)如cost_j a * sqrt(x_j) b或者擁堵成本凸函數(shù)。約束滿足所有區(qū)域的需求容量總和有上限投資總預(yù)算限制等。建模步驟定義集合候選站址集合J需求區(qū)域集合I。定義變量01變量y_j連續(xù)變量x_j以及可能的需求分配變量z_ij從站j滿足區(qū)域i的需求量。目標(biāo)函數(shù)最小化總成本 總建設(shè)成本非線性含y_j和x_j 總運(yùn)營/輸電成本可能是z_ij的函數(shù)。約束需求滿足對每個區(qū)域i∑_j z_ij demand_i。容量限制對每個站j∑_i z_ij x_j且x_j M * y_jBig-M法如果y_j0則x_j0。邏輯約束例如某個區(qū)域必須被至少一個站覆蓋∑_j a_ij * y_j 1其中a_ij表示站j是否能覆蓋區(qū)域i。預(yù)算約束總建設(shè)成本∑_j (F_j * y_j f(x_j)) ≤ Budget。求解策略這是一個混合整數(shù)非線性規(guī)劃問題。如果非線性部分可以線性化或分段線性化可以嘗試用Lingo或Matlab的intlinprog結(jié)合線性化技巧。如果非線性部分復(fù)雜可以考慮用啟發(fā)式算法如遺傳算法同時優(yōu)化y和x或者在y固定的情況下x的子問題是一個連續(xù)非線性規(guī)劃可以交替優(yōu)化。4.2 常見錯誤與調(diào)試心得在實(shí)現(xiàn)和求解這類模型時新手常會踩一些坑“無可行解”錯誤這是最令人頭疼的。首先檢查所有約束是否自相矛盾。例如邊界lb ub或者兩個約束聯(lián)合起來使得可行域?yàn)榭?。調(diào)試方法逐步注釋掉部分約束特別是非線性約束和復(fù)雜的邏輯約束先讓模型有解再逐個加入約束定位問題源。在Lingo中可以使用LINGO - Generate - Display model查看完整的線性化后的模型檢查約束。在Matlab中檢查A,b,Aeq,beq,lb,ub的維度是否正確?！敖獾馁|(zhì)量差”或“陷入局部最優(yōu)”對于非線性規(guī)劃這通常與初始點(diǎn)x0有關(guān)。策略進(jìn)行多初始點(diǎn)搜索。寫一個循環(huán)隨機(jī)生成多個x0確保在邊界內(nèi)或滿足簡單約束分別調(diào)用fmincon記錄最優(yōu)解。對于01規(guī)劃如果使用啟發(fā)式算法可以增加種群大小和迭代次數(shù)。Lingo報錯“NLP Solver failed”或求解時間過長非線性問題可能非凸Lingo的默認(rèn)本地求解器卡住了。嘗試在Lingo菜單LINGO - Options - Global Solver中勾選Use Global Solver。這會啟用全局優(yōu)化但計(jì)算時間會大幅增加。對于大規(guī)模問題這可能不現(xiàn)實(shí)此時需要考慮問題重構(gòu)或采用啟發(fā)式方法。數(shù)值不穩(wěn)定目標(biāo)函數(shù)或約束條件尺度差異巨大例如一個變量范圍是[0, 1]另一個是[10000, 20000]會導(dǎo)致求解器數(shù)值計(jì)算困難收斂緩慢或不準(zhǔn)確。解決方案對變量進(jìn)行縮放使其處于相近的數(shù)量級上例如都縮放到[0, 1]或[-1, 1]附近。這在Matlab和Lingo中都同樣重要。模型正確但求解器不收斂檢查優(yōu)化選項(xiàng)。在Matlab中適當(dāng)增加MaxIterations和MaxFunctionEvaluations。檢查收斂容差OptimalityTolerance和StepTolerance如果設(shè)置得太小可能永遠(yuǎn)達(dá)不到。有時候稍微放松容差如從1e-10調(diào)到1e-6就能讓求解器成功終止并獲得一個可接受的解。最后分享一個個人在數(shù)模競賽中處理復(fù)雜規(guī)劃問題的習(xí)慣永遠(yuǎn)從最簡單、最核心的模型版本開始。先忽略次要的非線性項(xiàng)用線性模型和01變量把主體邏輯跑通得到基準(zhǔn)解和運(yùn)行時間。然后再逐步加入非線性部分、更復(fù)雜的約束并觀察解的變化和計(jì)算時間的增長。這樣既能保證在有限時間內(nèi)有一個保底的模型也能有條理地評估模型復(fù)雜化帶來的收益與成本。記住在數(shù)模比賽中一個能跑出合理結(jié)果、邏輯清晰的簡化模型遠(yuǎn)勝過一個理論上完美但無法求解或求解不穩(wěn)定的復(fù)雜模型。