性規(guī)劃建模實(shí)戰(zhàn):從fmincon算法到投資組合優(yōu)化)
1. 項(xiàng)目概述為什么非線(xiàn)性規(guī)劃是建模的“硬骨頭”搞數(shù)學(xué)建模的朋友尤其是參加過(guò)國(guó)賽、美賽的應(yīng)該都深有體會(huì)線(xiàn)性規(guī)劃模型雖然基礎(chǔ)但真正讓你頭疼、讓你熬夜掉頭發(fā)的往往是那些“非線(xiàn)性”的家伙。題目里一旦出現(xiàn)“成本隨產(chǎn)量增加而邊際遞減”、“傳播速率與接觸人數(shù)成二次關(guān)系”、“收益與風(fēng)險(xiǎn)的非線(xiàn)性權(quán)衡”你的模型復(fù)雜度就立刻上了一個(gè)臺(tái)階。這就是非線(xiàn)性規(guī)劃Nonlinear Programming, NLP的領(lǐng)域。它不像線(xiàn)性規(guī)劃那樣有單純形法這種“萬(wàn)能鑰匙”其求解過(guò)程更像是在崎嶇的山地中尋找最高點(diǎn)你可能會(huì)陷入局部最優(yōu)的“小水坑”里而錯(cuò)過(guò)了遠(yuǎn)處的“山峰”。我最初接觸非線(xiàn)性規(guī)劃時(shí)也被各種算法、各種MATLAB函數(shù)搞得暈頭轉(zhuǎn)向。fmincon、fminunc、ga...每個(gè)函數(shù)都有一堆參數(shù)每個(gè)算法都有其適用場(chǎng)景用錯(cuò)了不僅結(jié)果不對(duì)還可能直接報(bào)錯(cuò)。市面上很多教程要么過(guò)于理論滿(mǎn)篇都是KKT條件、拉格朗日乘子要么過(guò)于簡(jiǎn)略只給個(gè)函數(shù)名了事。這篇內(nèi)容我就結(jié)合自己多年踩坑和輔導(dǎo)學(xué)生的經(jīng)驗(yàn)試圖把這塊“硬骨頭”啃碎、講透。我們的目標(biāo)很明確讓你在拿到一個(gè)非線(xiàn)性問(wèn)題后能迅速判斷其類(lèi)型選擇合適的MATLAB工具正確設(shè)置參數(shù)并解讀結(jié)果最終完成一篇邏輯清晰、求解可靠的建模論文。這篇文章將不僅僅是一個(gè)函數(shù)說(shuō)明書(shū)。我會(huì)從“道”與“術(shù)”兩個(gè)層面展開(kāi)“道”是理解非線(xiàn)性規(guī)劃問(wèn)題的本質(zhì)、分類(lèi)和求解思路“術(shù)”是掌握MATLAB特別是fmincon這個(gè)核心求解器的實(shí)戰(zhàn)技巧包括如何處理各種約束、如何選擇算法、如何調(diào)試和驗(yàn)證結(jié)果。無(wú)論你是建模新手還是想深化理解的進(jìn)階者希望都能從這里獲得可以直接“抄作業(yè)”的實(shí)用指南。2. 非線(xiàn)性規(guī)劃的核心思想與模型分類(lèi)在深入代碼之前我們必須把腦子里的概念理清楚。非線(xiàn)性規(guī)劃顧名思義就是目標(biāo)函數(shù)或約束條件中至少有一個(gè)是非線(xiàn)性函數(shù)的數(shù)學(xué)規(guī)劃問(wèn)題。它的通用形式可以寫(xiě)成最小化f(x)滿(mǎn)足A*x ≤ b(線(xiàn)性不等式約束)Aeq*x beq(線(xiàn)性等式約束)c(x) ≤ 0(非線(xiàn)性不等式約束)ceq(x) 0(非線(xiàn)性等式約束)lb ≤ x ≤ ub(決策變量上下界)其中x是決策變量向量f(x)是目標(biāo)函數(shù)c(x)和ceq(x)是非線(xiàn)性約束函數(shù)。注意A*x ≤ b和Aeq*x beq屬于線(xiàn)性約束它們可以和上下界lb,ub一起被高效處理。2.1 從幾何直觀理解“非線(xiàn)性”為什么非線(xiàn)性問(wèn)題難我們想象一個(gè)地形圖。線(xiàn)性規(guī)劃尋找最優(yōu)解就像在一個(gè)平坦的、傾斜的平面上找最低點(diǎn)你沿著坡度最陡的方向負(fù)梯度走一定能走到邊界上的最低點(diǎn)。而非線(xiàn)性規(guī)劃的地形可能是連綿起伏的山丘和山谷。凸函數(shù)與凹函數(shù)這是最重要的概念之一。如果目標(biāo)函數(shù)是凸函數(shù)并且可行域是凸集那么任何局部最優(yōu)解就是全局最優(yōu)解。這就像在一個(gè)碗狀地形里碗底只有一個(gè)。fmincon等求解器最喜歡這類(lèi)問(wèn)題求解效率高且結(jié)果可靠。反之如果是非凸的地形就像瑞士奶酪有很多局部最低點(diǎn)坑算法很容易掉進(jìn)某個(gè)坑里就停下來(lái)了以為找到了最優(yōu)解。約束的非線(xiàn)性約束條件畫(huà)出的可行域邊界不再是直線(xiàn)或平面可能是曲線(xiàn)或曲面。這會(huì)讓可行域的形狀變得復(fù)雜甚至不連通進(jìn)一步增加了搜索難度。2.2 常見(jiàn)非線(xiàn)性模型類(lèi)型舉例理解分類(lèi)才能對(duì)癥下藥。無(wú)約束非線(xiàn)性?xún)?yōu)化最簡(jiǎn)單的一類(lèi)只有目標(biāo)函數(shù)f(x)沒(méi)有約束。例如擬合一個(gè)非線(xiàn)性模型時(shí)的參數(shù)估計(jì)最小二乘法本質(zhì)就是在最小化誤差平方和這個(gè)非線(xiàn)性函數(shù)。MATLAB中用fminunc或fminsearch求解。僅含邊界約束變量有上下限比如物理意義要求數(shù)量非負(fù)x ≥ 0。這可以通過(guò)lb和ub參數(shù)輕松設(shè)置。線(xiàn)性約束非線(xiàn)性?xún)?yōu)化目標(biāo)函數(shù)非線(xiàn)性但所有約束都是線(xiàn)性的。這是fmincon非常擅長(zhǎng)處理的類(lèi)型。非線(xiàn)性約束優(yōu)化約束條件中出現(xiàn)了非線(xiàn)性函數(shù)。這是最復(fù)雜的一類(lèi)求解難度和計(jì)算量最大。例如在結(jié)構(gòu)設(shè)計(jì)中應(yīng)力、形變等約束常常是非線(xiàn)性的。二次規(guī)劃目標(biāo)函數(shù)是二次函數(shù)約束是線(xiàn)性的。它是非線(xiàn)性規(guī)劃中的一個(gè)特例有更高效的專(zhuān)門(mén)算法如quadprog。非線(xiàn)性最小二乘目標(biāo)函數(shù)形如一系列平方和的最小化。MATLAB提供了lsqnonlin和lsqcurvefit等專(zhuān)用函數(shù)比通用的fmincon效率更高。注意在數(shù)學(xué)建模競(jìng)賽中你遇到的大部分問(wèn)題都可以歸結(jié)為前三種。純非線(xiàn)性約束的問(wèn)題較少因?yàn)槠淝蠼夂捅硎鰧?duì)本科生而言挑戰(zhàn)較大。通常我們會(huì)嘗試通過(guò)變量代換、分段線(xiàn)性化等方法將非線(xiàn)性約束轉(zhuǎn)化為線(xiàn)性或邊界約束。2.3 求解的基本思路迭代與搜索所有數(shù)值優(yōu)化算法都是迭代法。它們從一個(gè)初始猜測(cè)解x0開(kāi)始然后根據(jù)某種規(guī)則算法產(chǎn)生一個(gè)搜索方向p_k和一個(gè)步長(zhǎng)α_k從而更新解x_{k1} x_k α_k * p_k。重復(fù)這個(gè)過(guò)程直到滿(mǎn)足停止條件如梯度足夠小、迭代次數(shù)超限、函數(shù)值變化不大等。不同的算法在于如何計(jì)算搜索方向p_k梯度下降法p_k取當(dāng)前點(diǎn)的負(fù)梯度方向。簡(jiǎn)單但可能在“山谷”中 zig-zag 前進(jìn)收斂慢。牛頓法利用目標(biāo)函數(shù)的二階導(dǎo)數(shù)Hessian矩陣信息能預(yù)測(cè)更優(yōu)的搜索方向收斂更快但計(jì)算 Hessian 矩陣代價(jià)高。擬牛頓法如BFGS牛頓法的改進(jìn)版通過(guò)迭代近似 Hessian 矩陣在收斂速度和計(jì)算成本間取得了很好的平衡。fmincon的內(nèi)點(diǎn)法和序列二次規(guī)劃算法都內(nèi)置了擬牛頓法更新。信賴(lài)域法在當(dāng)前位置的一個(gè)小“信賴(lài)域”內(nèi)用一個(gè)簡(jiǎn)單模型如二次模型近似原函數(shù)先優(yōu)化這個(gè)近似模型再根據(jù)結(jié)果調(diào)整信賴(lài)域大小和下一步迭代點(diǎn)。智能優(yōu)化算法如遺傳算法ga適用于高度非線(xiàn)性、非凸、甚至不連續(xù)的問(wèn)題。它們通過(guò)模擬自然進(jìn)化等過(guò)程進(jìn)行全局搜索不容易陷入局部最優(yōu)但計(jì)算量大且結(jié)果具有隨機(jī)性。對(duì)于建模我們的策略通常是先用智能算法如ga進(jìn)行全局粗略搜索找到潛力區(qū)域再將其結(jié)果作為fmincon的初始值x0進(jìn)行局部精細(xì)優(yōu)化。這能有效結(jié)合兩者的優(yōu)勢(shì)。3. MATLAB 核心求解器 fmincon 深度解析fmincon是MATLAB優(yōu)化工具箱中求解約束非線(xiàn)性多變量函數(shù)最小值的核心函數(shù)??梢哉f(shuō)掌握了fmincon就解決了80%以上的建模中的非線(xiàn)性規(guī)劃問(wèn)題。它的基本調(diào)用語(yǔ)法是[x, fval, exitflag, output, lambda, grad, hessian] fmincon(fun, x0, A, b, Aeq, beq, lb, ub, nonlcon, options)參數(shù)雖多但理解了邏輯就很簡(jiǎn)單。3.1 參數(shù)詳解與實(shí)戰(zhàn)準(zhǔn)備我們以一個(gè)經(jīng)典例子貫穿講解投資組合優(yōu)化。假設(shè)你有兩種資產(chǎn)期望收益率分別為r [0.15; 0.1]風(fēng)險(xiǎn)方差協(xié)方差矩陣為Sigma [0.2^2, 0.1*0.2*0.15; 0.1*0.2*0.15, 0.15^2]。你希望最小化投資組合的風(fēng)險(xiǎn)方差同時(shí)要求期望收益率不低于0.12且資金全部投出權(quán)重和為1每種資產(chǎn)投資比例非負(fù)。fun目標(biāo)函數(shù)句柄。這是需要你最小化的函數(shù)。它必須接受一個(gè)向量x并返回一個(gè)標(biāo)量。% 投資組合方差x * Sigma * x Sigma [0.04, 0.003; 0.003, 0.0225]; fun (x) x * Sigma * x;x0初始點(diǎn)。這是迭代的起點(diǎn)極其重要。一個(gè)糟糕的初始點(diǎn)可能導(dǎo)致算法收斂到局部最優(yōu)甚至失敗。對(duì)于有約束問(wèn)題x0必須是一個(gè)可行解即滿(mǎn)足所有約束。對(duì)于我們的例子可以設(shè)x0 [0.5; 0.5]它滿(mǎn)足權(quán)重和為1且非負(fù)。A, b線(xiàn)性不等式約束。表示A*x ≤ b。我們的收益率要求r*x ≥ 0.12是不等式需要轉(zhuǎn)化為-r*x ≤ -0.12。r [0.15; 0.1]; A -r; % 注意轉(zhuǎn)置因?yàn)?x 是列向量 b -0.12;Aeq, beq線(xiàn)性等式約束。表示Aeq*x beq。資金全部投出[1, 1] * x 1。Aeq [1, 1]; beq 1;lb, ub決策變量下界和上界。投資比例非負(fù)lb [0; 0]。沒(méi)有上限可以設(shè)為空[]。lb [0; 0]; ub []; % 表示無(wú)上界nonlcon非線(xiàn)性約束函數(shù)句柄。如果問(wèn)題有非線(xiàn)性約束c(x) ≤ 0或ceq(x) 0就需要定義這個(gè)函數(shù)。它接受x返回兩個(gè)向量[c, ceq]。本例沒(méi)有非線(xiàn)性約束設(shè)為[]。options優(yōu)化選項(xiàng)。這是高級(jí)用法和調(diào)試的關(guān)鍵。通過(guò)optimoptions(fmincon)創(chuàng)建。options optimoptions(fmincon, Display, iter, Algorithm, interior-point);Display, iter顯示每次迭代的詳細(xì)信息便于調(diào)試。Algorithm選擇核心算法這是重中之重下文詳述。3.2 算法選擇四大內(nèi)功心法fmincon提供了多種算法對(duì)應(yīng)不同的“內(nèi)功心法”。選擇不當(dāng)輕則效率低下重則無(wú)法收斂。interior-point內(nèi)點(diǎn)法默認(rèn)算法原理從可行域內(nèi)部出發(fā)通過(guò)構(gòu)造障礙函數(shù)在迭代過(guò)程中始終保持在可行域內(nèi)部并逐漸逼近邊界上的最優(yōu)解。優(yōu)點(diǎn)處理大規(guī)模問(wèn)題變量多、約束多性能優(yōu)秀特別擅長(zhǎng)處理不等式約束和邊界約束。對(duì)于我們的投資組合問(wèn)題有不等式和邊界約束它是很好的選擇。缺點(diǎn)對(duì)于問(wèn)題尺度較小或主要包含等式約束的問(wèn)題可能不是最快。適用默認(rèn)首選尤其當(dāng)你的問(wèn)題包含大量不等式約束時(shí)。sqp序列二次規(guī)劃法原理在每一步迭代用二次函數(shù)近似目標(biāo)函數(shù)用線(xiàn)性函數(shù)近似約束求解一個(gè)二次規(guī)劃子問(wèn)題從而確定搜索方向。優(yōu)點(diǎn)對(duì)于中小規(guī)模問(wèn)題特別是非線(xiàn)性約束問(wèn)題往往非常高效和精確。它能很好地利用目標(biāo)函數(shù)和約束的函數(shù)值、梯度信息。缺點(diǎn)對(duì)于大規(guī)模問(wèn)題子問(wèn)題的求解可能變得昂貴。適用問(wèn)題規(guī)模不大變量數(shù)幾百以?xún)?nèi)且含有非線(xiàn)性約束時(shí)可以?xún)?yōu)先嘗試sqp。active-set有效集法原理猜測(cè)哪些約束在最優(yōu)解處是“激活”的等式成立然后主要在這些約束構(gòu)成的子空間上進(jìn)行優(yōu)化。優(yōu)點(diǎn)能提供非常精確的拉格朗日乘子估計(jì)lambda輸出對(duì)于需要靈敏度分析的情況很有用。缺點(diǎn)不適合大規(guī)模問(wèn)題迭代過(guò)程中可能需要在有效集之間頻繁切換效率可能不如內(nèi)點(diǎn)法。適用需要高精度的乘子信息或者問(wèn)題規(guī)模較小且已知最優(yōu)解大概在哪些約束邊界上時(shí)。trust-region-reflective信賴(lài)域反射法原理屬于信賴(lài)域法要求目標(biāo)函數(shù)能提供梯度并且約束只能是邊界約束或線(xiàn)性等式約束。不能處理非線(xiàn)性約束或線(xiàn)性不等式約束。優(yōu)點(diǎn)對(duì)于邊界約束或線(xiàn)性等式約束的問(wèn)題如果提供了梯度此法可能非常高效。缺點(diǎn)適用范圍窄。適用只有邊界約束或邊界約束線(xiàn)性等式約束的問(wèn)題且你能計(jì)算或提供目標(biāo)函數(shù)的梯度。實(shí)操心得對(duì)于建模競(jìng)賽中的大部分問(wèn)題我的建議是無(wú)腦先用interior-point。如果求解失敗或結(jié)果可疑再?lài)L試sqp。除非問(wèn)題有特殊結(jié)構(gòu)如純邊界約束否則很少需要手動(dòng)切到其他算法。將Display設(shè)為iter觀察迭代過(guò)程是判斷算法是否正常工作的好方法。3.3 完整求解示例與結(jié)果解讀現(xiàn)在我們把所有部分組合起來(lái)求解投資組合問(wèn)題。% 1. 定義問(wèn)題數(shù)據(jù) Sigma [0.04, 0.003; 0.003, 0.0225]; % 協(xié)方差矩陣 r [0.15; 0.1]; % 期望收益率 targetReturn 0.12; % 目標(biāo)最低收益率 % 2. 定義目標(biāo)函數(shù) fun (x) x * Sigma * x; % 3. 初始點(diǎn) (一個(gè)可行的猜測(cè)) x0 [0.5; 0.5]; % 4. 線(xiàn)性不等式約束: r*x targetReturn - -r*x -targetReturn A -r; b -targetReturn; % 5. 線(xiàn)性等式約束: 權(quán)重之和為1 Aeq [1, 1]; beq 1; % 6. 邊界約束: 權(quán)重非負(fù) lb [0; 0]; ub []; % 無(wú)上界 % 7. 非線(xiàn)性約束: 無(wú) nonlcon []; % 8. 設(shè)置選項(xiàng)使用內(nèi)點(diǎn)法并顯示迭代信息 options optimoptions(fmincon, Display, iter, Algorithm, interior-point); % 9. 調(diào)用 fmincon 求解 [x_opt, fval_opt, exitflag, output, lambda] fmincon(fun, x0, A, b, Aeq, beq, lb, ub, nonlcon, options); % 10. 輸出結(jié)果 fprintf(最優(yōu)投資比例\n); fprintf( 資產(chǎn)1: %.4f\n, x_opt(1)); fprintf( 資產(chǎn)2: %.4f\n, x_opt(2)); fprintf(投資組合最小方差風(fēng)險(xiǎn): %.6f\n, fval_opt); fprintf(投資組合期望收益率: %.4f\n, r*x_opt); fprintf(退出標(biāo)志 exitflag: %d\n, exitflag); fprintf(迭代次數(shù): %d\n, output.iterations); fprintf(函數(shù)計(jì)算次數(shù): %d\n, output.funcCount);運(yùn)行后你會(huì)在命令窗口看到詳細(xì)的迭代過(guò)程最后得到結(jié)果。關(guān)鍵是要會(huì)看exitflagexitflag 0優(yōu)化成功收斂。通常是1表示一階最優(yōu)性條件在指定容差內(nèi)滿(mǎn)足。exitflag 0迭代次數(shù)或函數(shù)計(jì)算次數(shù)超過(guò)了options.MaxIterations或options.MaxFunctionEvaluations。此時(shí)結(jié)果可能不是最優(yōu)的需要增加迭代上限或檢查問(wèn)題。exitflag 0優(yōu)化失敗。常見(jiàn)的有 -2未找到可行點(diǎn)檢查初始點(diǎn)x0和約束-1被輸出函數(shù)或繪圖函數(shù)終止。在我們的例子中應(yīng)該會(huì)得到exitflag 1以及類(lèi)似x_opt [0.4; 0.6]的結(jié)果。lambda結(jié)構(gòu)體包含了約束對(duì)應(yīng)的拉格朗日乘子其eqlin字段對(duì)應(yīng)等式約束的乘子ineqlin對(duì)應(yīng)不等式約束lower/upper對(duì)應(yīng)邊界約束。乘子的絕對(duì)值大小反映了該約束的“緊度”或“價(jià)值”。4. 從理論到實(shí)戰(zhàn)復(fù)雜模型構(gòu)建與求解策略掌握了基礎(chǔ)我們來(lái)看幾個(gè)建模中更典型的復(fù)雜場(chǎng)景及其處理策略。4.1 場(chǎng)景一目標(biāo)函數(shù)或約束需要外部數(shù)據(jù)或復(fù)雜計(jì)算很多時(shí)候你的目標(biāo)函數(shù)f(x)不是一個(gè)簡(jiǎn)單的數(shù)學(xué)表達(dá)式而是一個(gè)“黑箱”過(guò)程。比如x是某個(gè)系統(tǒng)的設(shè)計(jì)參數(shù)f(x)是需要調(diào)用一個(gè)仿真程序如 Simulink 模型、有限元分析才能計(jì)算出的性能指標(biāo)。策略封裝函數(shù)將仿真過(guò)程寫(xiě)成一個(gè)獨(dú)立的MATLAB函數(shù)mySimulation(x)該函數(shù)接受參數(shù)x運(yùn)行仿真并返回標(biāo)量結(jié)果如最大應(yīng)力、總能耗。然后fun句柄指向這個(gè)封裝函數(shù)。function cost myComplexObjective(x) % x 是設(shè)計(jì)參數(shù) % 1. 根據(jù)x設(shè)置模型參數(shù) setModelParameters(x); % 2. 運(yùn)行外部仿真或復(fù)雜計(jì)算這里用耗時(shí)計(jì)算模擬 result runExternalSimulation(); % 假設(shè)這個(gè)函數(shù)很耗時(shí) % 3. 從結(jié)果中提取目標(biāo)值 cost extractCostFromResult(result); end % 在優(yōu)化中調(diào)用 fun myComplexObjective;注意事項(xiàng)這類(lèi)問(wèn)題計(jì)算一次目標(biāo)函數(shù)代價(jià)很高。務(wù)必設(shè)置合理的options.MaxFunctionEvaluations和options.MaxIterations避免無(wú)意義的長(zhǎng)時(shí)間運(yùn)行。同時(shí)考慮使用UseParallel選項(xiàng)為true如果目標(biāo)函數(shù)計(jì)算可以并行化的話(huà)能極大加速。4.2 場(chǎng)景二多目標(biāo)優(yōu)化問(wèn)題現(xiàn)實(shí)中我們常需要權(quán)衡多個(gè)目標(biāo)例如“成本最低”和“性能最好”。這被稱(chēng)為多目標(biāo)優(yōu)化其解不是一個(gè)點(diǎn)而是一個(gè)“帕累托前沿”P(pán)areto Front。策略加權(quán)求和法最直接的方法是將多目標(biāo)轉(zhuǎn)化為單目標(biāo)F(x) w1 * f1(x) w2 * f2(x)。通過(guò)調(diào)整權(quán)重w1,w2可以得到前沿上的不同點(diǎn)。w1 0.7; w2 0.3; % 權(quán)重代表決策者的偏好 fun (x) w1 * costFunction(x) w2 * (-performanceFunction(x)); % 注意性能可能是最大化加負(fù)號(hào)轉(zhuǎn)為最小化策略目標(biāo)規(guī)劃法設(shè)定一個(gè)理想的目標(biāo)值然后最小化與它的偏差。targetCost 1000; targetPerf 50; fun (x) abs(costFunction(x) - targetCost) abs(performanceFunction(x) - targetPerf);策略使用帕累托搜索算法對(duì)于復(fù)雜的多目標(biāo)問(wèn)題可以使用MATLAB的paretosearch或gamultiobj多目標(biāo)遺傳算法來(lái)直接尋找近似帕累托前沿。這在建模中是非常高級(jí)和出彩的技巧。4.3 場(chǎng)景三含非線(xiàn)性約束的問(wèn)題假設(shè)在投資組合中我們?cè)黾右粋€(gè)非線(xiàn)性約束要求兩種資產(chǎn)權(quán)重的乘積不超過(guò)某個(gè)值模擬某種關(guān)聯(lián)性限制即x1 * x2 ≤ 0.2。這時(shí)就需要定義nonlcon函數(shù)。function [c, ceq] myNonlcon(x) % 非線(xiàn)性不等式約束 c(x) 0 c x(1) * x(2) - 0.2; % 注意要求 c 0所以是 x1*x2 - 0.2 0 % 非線(xiàn)性等式約束 ceq(x) 0 ceq []; % 本例沒(méi)有非線(xiàn)性等式約束 end % 在 fmincon 調(diào)用中傳入 nonlcon myNonlcon; [x_opt, fval] fmincon(fun, x0, A, b, Aeq, beq, lb, ub, myNonlcon, options);踩坑記錄非線(xiàn)性約束函數(shù)的編寫(xiě)是錯(cuò)誤高發(fā)區(qū)。務(wù)必記住c(x) ≤ 0和ceq(x) 0。經(jīng)常有人把不等式方向?qū)懛?。另外確保nonlcon函數(shù)能正確返回兩個(gè)輸出[c, ceq]即使其中一個(gè)為空。4.4 場(chǎng)景四變量離散或整數(shù)規(guī)劃如果變量只能取整數(shù)如設(shè)備臺(tái)數(shù)或離散值如標(biāo)準(zhǔn)尺寸問(wèn)題就變成了混合整數(shù)非線(xiàn)性規(guī)劃。fmincon無(wú)法直接處理。策略連續(xù)松弛圓整先忽略整數(shù)約束用fmincon求解連續(xù)問(wèn)題。得到連續(xù)最優(yōu)解后將其圓整到最近的整數(shù)或離散值。但要注意圓整后的解可能不可行違反約束或遠(yuǎn)離真正的最優(yōu)解。這只是一種近似啟發(fā)式方法。策略使用專(zhuān)用求解器MATLAB的全局優(yōu)化工具箱提供了ga遺傳算法支持整數(shù)約束和surrogateopt代理優(yōu)化等可以直接處理整數(shù)變量。對(duì)于復(fù)雜的整數(shù)非線(xiàn)性規(guī)劃可能需要更專(zhuān)業(yè)的工具如intlinprog僅線(xiàn)性或第三方求解器。5. 調(diào)試、驗(yàn)證與結(jié)果分析避免“垃圾進(jìn)垃圾出”優(yōu)化求解器不是魔法它只是忠實(shí)地執(zhí)行你定義的模型。如果模型有誤、初始點(diǎn)太差或參數(shù)設(shè)置不當(dāng)?shù)玫降慕Y(jié)果就是無(wú)意義的。因此求解后的調(diào)試和驗(yàn)證至關(guān)重要。5.1 診斷求解失敗如果exitflag不是正數(shù)按以下步驟排查檢查初始點(diǎn)x0它必須是可行的用x0代入所有約束條件驗(yàn)算。對(duì)于不等式A*x0 b和c(x0) 0以及等式Aeq*x0 beq和ceq(x0) 0在容差內(nèi)。一個(gè)簡(jiǎn)單的方法是先求解一個(gè)可行性問(wèn)題或者手動(dòng)構(gòu)造一個(gè)明顯的可行點(diǎn)。檢查約束矛盾約束是否可能相互沖突導(dǎo)致無(wú)解例如兩個(gè)不等式約束可能定義了空集??梢試L試放松或移除部分約束看問(wèn)題是否變得可行。檢查梯度/導(dǎo)數(shù)信息如果你通過(guò)options提供了梯度或 Hessian 函數(shù)SpecifyObjectiveGradient,true務(wù)必檢查其計(jì)算是否正確。一個(gè)錯(cuò)誤的梯度會(huì)導(dǎo)致算法在錯(cuò)誤的方向搜索??梢允褂胏heckGradients選項(xiàng)或fmincon的CheckGradients選項(xiàng)進(jìn)行數(shù)值驗(yàn)證。調(diào)整算法和選項(xiàng)換一個(gè)算法試試如從interior-point換到sqp。增加迭代次數(shù)和函數(shù)計(jì)算次數(shù)上限MaxIterations,MaxFunctionEvaluations。放寬最優(yōu)性容差OptimalityTolerance或約束容差ConstraintTolerance尤其是在目標(biāo)函數(shù)或約束值非常小或非常大時(shí)。嘗試不同的初始點(diǎn)x0。多跑幾次從隨機(jī)初始點(diǎn)開(kāi)始觀察是否收斂到同一點(diǎn)。5.2 驗(yàn)證最優(yōu)解即使exitflag 0也未必是全局最優(yōu)尤其是對(duì)于非凸問(wèn)題??尚行则?yàn)證將最優(yōu)解x_opt代回所有約束確保滿(mǎn)足在ConstraintTolerance內(nèi)。局部最優(yōu)性檢查觀察output.firstorderopt輸出它是一階最優(yōu)性條件的度量值越小越好接近OptimalityTolerance。對(duì)于無(wú)約束問(wèn)題可以手動(dòng)計(jì)算梯度gradient(fun, x_opt)看其范數(shù)是否接近零。敏感性分析拉格朗日乘子lambda結(jié)構(gòu)體中的乘子提供了寶貴信息。對(duì)于一個(gè)不等式約束如果其乘子lambda.ineqlin(i)的絕對(duì)值很大說(shuō)明這個(gè)約束是“緊”的活躍的放松它會(huì)對(duì)目標(biāo)函數(shù)有顯著改善。如果乘子為0則該約束在最優(yōu)解處不活躍。全局最優(yōu)性試探多初始點(diǎn)法從多個(gè)隨機(jī)初始點(diǎn)運(yùn)行fmincon比較得到的目標(biāo)函數(shù)值。如果都收斂到相同或相近的值則全局最優(yōu)的可能性增大。使用全局優(yōu)化求解器用ga遺傳算法或particleswarm粒子群算法等全局優(yōu)化器求解同一個(gè)問(wèn)題。比較它們找到的最佳值與fmincon的結(jié)果。如果fmincon的結(jié)果差很多說(shuō)明它可能陷入了局部最優(yōu)。此時(shí)可以將ga找到的解作為fmincon的初始點(diǎn)進(jìn)行“雜交”優(yōu)化。5.3 結(jié)果呈現(xiàn)與論文寫(xiě)作在建模論文中不能只扔出一個(gè)數(shù)字。清晰表述模型用數(shù)學(xué)公式明確寫(xiě)出目標(biāo)函數(shù)和所有約束。說(shuō)明求解工具寫(xiě)明“使用MATLAB R2023a的優(yōu)化工具箱中的fmincon函數(shù)進(jìn)行求解采用內(nèi)點(diǎn)算法”。報(bào)告關(guān)鍵參數(shù)給出初始點(diǎn)x0、重要的options設(shè)置如算法、容差。展示求解結(jié)果以表格形式呈現(xiàn)最優(yōu)解x_opt、最優(yōu)目標(biāo)值fval、關(guān)鍵約束的滿(mǎn)足情況。進(jìn)行分析討論靈敏度分析改變模型中的某個(gè)參數(shù)如投資組合中的目標(biāo)收益率targetReturn重新求解觀察最優(yōu)解如何變化??梢岳L制出“有效前沿”曲線(xiàn)風(fēng)險(xiǎn) vs 收益。模型穩(wěn)健性如果數(shù)據(jù)有微小波動(dòng)最優(yōu)解變化大嗎可以通過(guò)在參數(shù)上加微小擾動(dòng)來(lái)測(cè)試。結(jié)果解釋最優(yōu)解在現(xiàn)實(shí)中有何意義權(quán)重分配是否符合直覺(jué)如果不符合是模型漏掉了什么重要約束嗎6. 高級(jí)技巧與性能優(yōu)化當(dāng)問(wèn)題規(guī)模變大或函數(shù)計(jì)算昂貴時(shí)這些技巧能幫你節(jié)省大量時(shí)間。6.1 提供解析梯度與Hessian默認(rèn)情況下fmincon使用有限差分法數(shù)值估算梯度。這需要多次調(diào)用目標(biāo)函數(shù)且精度有限。如果你能提供目標(biāo)函數(shù)梯度的解析表達(dá)式能極大提升速度和精度。function [f, g] myObjectiveWithGradient(x) % 計(jì)算目標(biāo)函數(shù)值 f f x(1)^2 sin(x(2)); % 計(jì)算梯度 g [df/dx1; df/dx2] if nargout 1 % 只有當(dāng)需要梯度時(shí)才計(jì)算 g [2*x(1); cos(x(2))]; end end options optimoptions(fmincon, SpecifyObjectiveGradient, true); [x, fval] fmincon(myObjectiveWithGradient, x0, ..., options);對(duì)于Hessian矩陣也是如此HessianFcn。對(duì)于大規(guī)模問(wèn)題提供梯度收益顯著。6.2 并行計(jì)算加速如果目標(biāo)函數(shù)或約束函數(shù)的計(jì)算可以向量化或獨(dú)立進(jìn)行開(kāi)啟并行計(jì)算能成倍縮短時(shí)間。options optimoptions(fmincon, UseParallel, true);在調(diào)用fmincon前確保已經(jīng)通過(guò)parpool開(kāi)啟了并行池。這特別適用于前述的“黑箱”仿真類(lèi)目標(biāo)函數(shù)或者使用多初始點(diǎn)法時(shí)。6.3 變量縮放與預(yù)處理優(yōu)化問(wèn)題的“條件數(shù)”很重要。如果變量x1的范圍是[0, 1]而x2的范圍是[1000, 2000]這會(huì)導(dǎo)致數(shù)值問(wèn)題使算法收斂緩慢。策略縮放變量。引入新的縮放變量y使得x scale * y讓y的各分量量級(jí)大致相當(dāng)。例如令y1 x1,y2 x2 / 1000。在目標(biāo)函數(shù)和約束中都用y來(lái)表示最后結(jié)果再轉(zhuǎn)換回x。這能顯著改善算法的數(shù)值穩(wěn)定性。6.4 利用問(wèn)題結(jié)構(gòu)稀疏性與對(duì)稱(chēng)性對(duì)于大規(guī)模問(wèn)題如果 Jacobian 矩陣約束的導(dǎo)數(shù)或 Hessian 矩陣是稀疏的一定要通過(guò)options告知求解器JacobPattern,HessPattern這能節(jié)省大量?jī)?nèi)存和計(jì)算時(shí)間。在建模競(jìng)賽的超大規(guī)模問(wèn)題中這一點(diǎn)可能至關(guān)重要。7. 常見(jiàn)問(wèn)題與排查技巧實(shí)錄這里匯總了我自己和學(xué)生們?cè)趯?shí)戰(zhàn)中踩過(guò)的坑和解決方法。問(wèn)題現(xiàn)象可能原因排查與解決思路exitflag -2(找不到可行點(diǎn))1. 初始點(diǎn)x0不可行。2. 約束條件相互矛盾可行域?yàn)榭铡?.驗(yàn)證x0將其代入所有約束計(jì)算。手動(dòng)構(gòu)造一個(gè)簡(jiǎn)單的可行點(diǎn)如所有邊界的中點(diǎn)。2.松弛約束暫時(shí)注釋掉部分約束特別是非線(xiàn)性約束看問(wèn)題是否變得可行。逐步添加約束以定位矛盾點(diǎn)。exitflag 0(達(dá)到迭代上限)1. 問(wèn)題太復(fù)雜需要更多迭代。2. 算法在平緩區(qū)域“蠕動(dòng)”收斂慢。1.增加限制options.MaxIterations和options.MaxFunctionEvaluations。2.檢查收斂趨勢(shì)設(shè)置Display, iter看目標(biāo)函數(shù)值是否還在穩(wěn)定下降。如果下降緩慢可能是接近最優(yōu)解可以適當(dāng)收緊OptimalityTolerance或StepTolerance以提前停止。3.更換算法或提供梯度。結(jié)果對(duì)初始點(diǎn)x0敏感問(wèn)題是非凸的存在多個(gè)局部最優(yōu)解。1.多初始點(diǎn)法用MultiStart或GlobalSearch封裝fmincon自動(dòng)從多個(gè)初始點(diǎn)搜索。2.使用全局優(yōu)化器先用ga進(jìn)行全局搜索再用其結(jié)果作為fmincon的初始點(diǎn)。求解速度極慢1. 目標(biāo)/約束函數(shù)計(jì)算耗時(shí)。2. 問(wèn)題規(guī)模大。3. 數(shù)值條件差變量尺度差異大。1.提供解析導(dǎo)數(shù)梯度、Hessian。2.開(kāi)啟并行計(jì)算UseParallel, true。3.進(jìn)行變量縮放。4. 嘗試更高效的算法如對(duì)邊界約束問(wèn)題用trust-region-reflective。得到的結(jié)果明顯不合理(如負(fù)的投資比例)1. 邊界約束lb設(shè)置錯(cuò)誤或未設(shè)置。2. 模型本身有誤如目標(biāo)函數(shù)符號(hào)反了。3. 算法陷入了一個(gè)很差的局部最優(yōu)。1.仔細(xì)檢查模型打印出目標(biāo)函數(shù)和約束在最優(yōu)解處的值手動(dòng)驗(yàn)算。2.檢查邊界確認(rèn)lb和ub是否正確施加。3.從不同初始點(diǎn)重新求解對(duì)比結(jié)果。4.簡(jiǎn)化問(wèn)題先去掉復(fù)雜約束看基礎(chǔ)版本是否合理。fmincon提示“用戶(hù)提供的目標(biāo)函數(shù)返回了NaN或Inf”目標(biāo)函數(shù)或約束函數(shù)在某些x處未定義如除以零、對(duì)負(fù)數(shù)取對(duì)數(shù)。1. 在函數(shù)內(nèi)部添加防御性代碼檢查輸入x的有效性對(duì)非法操作返回一個(gè)很大的懲罰值如1e10引導(dǎo)優(yōu)化器離開(kāi)該區(qū)域。2. 收緊變量的上下界lb,ub避免函數(shù)未定義的區(qū)域。最后再分享一個(gè)我常用的調(diào)試流程從簡(jiǎn)到繁逐步驗(yàn)證。先構(gòu)建一個(gè)最簡(jiǎn)單的、有已知解析解或明顯答案的模型用fmincon求解確認(rèn)代碼框架和模型表述正確。然后逐步添加復(fù)雜的約束和非線(xiàn)性項(xiàng)每加一步都驗(yàn)證結(jié)果的合理性。這樣能最快地定位問(wèn)題所在避免在復(fù)雜的模型里大海撈針。非線(xiàn)性規(guī)劃求解就像偵探破案需要邏輯、耐心和對(duì)細(xì)節(jié)的把握。希望這篇超詳細(xì)的指南能成為你建模工具箱里一件稱(chēng)手的利器。