數(shù)學(xué)建模-非線性規(guī)劃模型用MATLAB++LINGO

約束非線性規(guī)劃非現(xiàn)性規(guī)劃的基本概念

定義如果目標(biāo)函數(shù)或約束條件中至少有一個是非線性函數(shù)時的最優(yōu)化問題就叫做非線性規(guī)劃問題．

一般形式:

（1）其中，是定義在En上的實值函數(shù)，簡記:

其它情況:

求目標(biāo)函數(shù)的最大值或約束條件為小于等于零的情況，都可通過取其相反數(shù)化為上述一般形式．

定義1把滿足問題（1）中條件的解稱為可行解（或可行點），所有可行點的集合稱為可行集（或可行域）．記為D．即問題(1)可簡記為．定義2

對于問題(1)，設(shè)，若存在,使得對一切，且，都有，則稱X*是f(X)在D上的局部極小值點（局部最優(yōu)解）．特別地當(dāng)時，若則稱X*是f(X)在D上的嚴(yán)格局部極小值點（嚴(yán)格局部最優(yōu)解）．定義3對于問題(1),設(shè),對任意的,都有則稱X*是f(X)在D上的全局極小值點（全局最優(yōu)解）．特別地當(dāng)時，若，則稱X*是f(X)在D上的嚴(yán)格全局極小值點（嚴(yán)格全局最優(yōu)解）．

返回非線性規(guī)劃的基本解法SUTM外點法SUTM內(nèi)點法（障礙罰函數(shù)法）1、罰函數(shù)法2、近似規(guī)劃法

返回

罰函數(shù)法

罰函數(shù)法基本思想是通過構(gòu)造罰函數(shù)把約束問題轉(zhuǎn)化為一系列無約束最優(yōu)化問題，進(jìn)而用無約束最優(yōu)化方法去求解．這類方法稱為序列無約束最小化方法(SequentialUnconstrainedMinizationTechnique)．簡稱為SUMT法．其一為SUMT外點法，其二為SUMT內(nèi)點法．

其中T(X,M)稱為罰函數(shù)，M稱為罰因子，帶M的項稱為罰項，這里的罰函數(shù)只對不滿足約束條件的點實行懲罰：當(dāng)時，滿足各，故罰項=0，不受懲罰．當(dāng)時，必有的約束條件，故罰項>0，要受懲罰．SUTM外點法

罰函數(shù)法的缺點是：每個近似最優(yōu)解Xk往往不是容許解，而只能近似滿足約束，在實際問題中這種結(jié)果可能不能使用；在解一系列無約束問題中，計算量太大，特別是隨著Mk的增大，可能導(dǎo)致錯誤．1、任意給定初始點X0，取M1>1，給定允許誤差，令k=1；2、求無約束極值問題的最優(yōu)解，設(shè)為Xk=X(Mk)，即；3、若存在，使，則取Mk>M()令k=k+1返回（2），否則，停止迭代．得最優(yōu)解.計算時也可將收斂性判別準(zhǔn)則改為.SUTM外點法(罰函數(shù)法)的迭代步驟例題.minx2

s.t.x>=1顯然本問題的最優(yōu)解為x*=1用SMT外點法:

T(x,M)=x2+M[min(0,x-1)]2

=求minT(x,M).本題可由

T’(x,M)=2x+2M(x-1)=0,解得:x=M/(1+M),M趨于無窮.可知x從小于1趨于1,罰函數(shù)從外部趨于最優(yōu)解.

SUTM內(nèi)點法（障礙函數(shù)法）

內(nèi)點法的迭代步驟例題.minx2

s.t.x>=1顯然本問題的最優(yōu)解為x*=1用SMT內(nèi)點法:罰函數(shù)I(x,rk)=x2

–rkln(x-1),求其最值.可見,x從大于1趨于1.

近似規(guī)劃法的基本思想：將問題(3)中的目標(biāo)函數(shù)和約束條件近似為線性函數(shù)，并對變量的取值范圍加以限制，從而得到一個近似線性規(guī)劃問題，再用單純形法求解之，把其符合原始條件的最優(yōu)解作為(3)的解的近似．近似規(guī)劃法每得到一個近似解后，都從這點出發(fā)，重復(fù)以上步驟．

這樣，通過求解一系列線性規(guī)劃問題，產(chǎn)生一個由線性規(guī)劃最優(yōu)解組成的序列，經(jīng)驗表明，這樣的序列往往收斂于非線性規(guī)劃問題的解。

近似規(guī)劃法的算法步驟如下

返回例題1、用近似規(guī)劃法求解下列問題。解：第一次迭代：在點x1處將g1(x),g2(x)線性化。步長限制：解（1）（2）（3）（4），得，檢驗，可得x2不滿足原始約束。第二次迭代，減少解（1,2,3,5)，得x2=(4,3)T,f(x2)=-11.所以x2為最小值點.用MATLAB軟件求解,其輸入格式如下:

1. x=quadprog(H,C,A,b);2. x=quadprog(H,C,A,b,Aeq,beq);3. x=quadprog(H,C,A,b,Aeq,beq,VLB,VUB);4. x=quadprog(H,C,A,b,Aeq,beq,VLB,VUB,X0);5. x=quadprog(H,C,A,b,Aeq,beq,VLB,VUB,X0,options);6. [x,fval]=quaprog(...);7. [x,fval,exitflag]=quaprog(...);8. [x,fval,exitflag,output]=quaprog(...);1、二次規(guī)劃例1

minf(x1,x2)=-2x1-6x2+x12-2x1x2+2x22

s.t.x1+x2≤2-x1+2x2≤2x1≥0,x2≥01、寫成標(biāo)準(zhǔn)形式：2、輸入命令：

H=[1-1;-12];c=[-2;-6];A=[11;-12];b=[2;2];

Aeq=[];beq=[];VLB=[0;0];VUB=[];[x,z]=quadprog(H,c,A,b,Aeq,beq,VLB,VUB)3、運(yùn)算結(jié)果為：

x=0.66671.3333z=-8.2222s.t.

1.首先建立M文件fun.m,定義目標(biāo)函數(shù)F（X）:functionf=fun(X);f=F(X);2、一般非線性規(guī)劃

其中X為n維變元向量，G(X)與Ceq(X)均為非線性函數(shù)組成的向量，其它變量的含義與線性規(guī)劃、二次規(guī)劃中相同.用Matlab求解上述問題，基本步驟分三步：3.建立主程序.非線性規(guī)劃求解的函數(shù)是fmincon,命令的基本格式如下：

(1)x=fmincon(‘fun’,X0,A,b)

(2)x=fmincon(‘fun’,X0,A,b,Aeq,beq)

(3)x=fmincon(‘fun’,X0,A,b,Aeq,beq,VLB,VUB)

(4)x=fmincon(‘fun’,X0,A,b,Aeq,beq,VLB,VUB,’nonlcon’)(5)x=fmincon(‘fun’,X0,A,b,Aeq,beq,VLB,VUB,’nonlcon’,options)

(6)[x,fval]=

fmincon(...)

(7)[x,fval,exitflag]=

fmincon(...)(8)[x,fval,exitflag,output]=

fmincon(...)輸出極值點M文件迭代的初值參數(shù)說明變量上下限注意：[1]fmincon函數(shù)提供了大型優(yōu)化算法和中型優(yōu)化算法。默認(rèn)時，若在fun函數(shù)中提供了梯度（options參數(shù)的GradObj設(shè)置為’on’），并且只有上下界存在或只有等式約束，fmincon函數(shù)將選擇大型算法。當(dāng)既有等式約束又有梯度約束時，使用中型算法。[2]fmincon函數(shù)的中型算法使用的是序列二次規(guī)劃法。在每一步迭代中求解二次規(guī)劃子問題，并用BFGS法更新拉格朗日Hessian矩陣。[3]fmincon函數(shù)可能會給出局部最優(yōu)解，這與初值X0的選取有關(guān)。1、寫成標(biāo)準(zhǔn)形式：

s.t.

2x1+3x26

s.tx1+4x25x1,x20例22、先建立M-文件fun3.m:

functionf=fun3(x);f=-x(1)-2*x(2)+(1/2)*x(1)^2+(1/2)*x(2)^23、再建立主程序youh2.m：

x0=[1;1];A=[23;14];b=[6;5];

Aeq=[];beq=[];VLB=[0;0];VUB=[];[x,fval]=fmincon('fun3',x0,A,b,Aeq,beq,VLB,VUB)4、運(yùn)算結(jié)果為：

x=0.76471.0588

fval=-2.02941．先建立M文件fun4.m,定義目標(biāo)函數(shù):

functionf=fun4(x);f=exp(x(1))*(4*x(1)^2+2*x(2)^2+4*x(1)*x(2)+2*x(2)+1);例32．再建立M文件mycon.m定義非線性約束：

function[g,ceq]=mycon(x)g=x(1)+x(2);

ceq=[1.5+x(1)*x(2)-x(1)-x(2);-x(1)*x(2)-10];3．主程序youh3.m為:x0=[-1;1];A=[];b=[];Aeq=[11];beq=[0];vlb=[];vub=[];[x,fval]=fmincon('fun4',x0,A,b,Aeq,beq,vlb,vub,'mycon')3.運(yùn)算結(jié)果為：

x=-1.22501.2250

fval=1.8951

例4

1．先建立M-文件fun.m定義目標(biāo)函數(shù):functionf=fun(x);f=-2*x(1)-x(2);2．再建立M文件mycon2.m定義非線性約束：

function[g,ceq]=mycon2(x)g=[x(1)^2+x(2)^2-25;x(1)^2-x(2)^2-7];

3.主程序fxx.m為:x0=[3;2.5];VLB=[00];VUB=[510];[x,fval,exitflag,output]=fmincon('fun',x0,[],[],[],[],VLB,VUB,'mycon2')4.運(yùn)算結(jié)果為:x=4.00003.0000fval=-11.0000exitflag=1output=iterations:4

funcCount:17

stepsize:1algorithm:[1x44char]

firstorderopt:[]

cgiterations:[]應(yīng)用實例：供應(yīng)與選址

某公司有6個建筑工地要開工，每個工地的位置（用平面坐標(biāo)系a，b表示，距離單位：千米）及水泥日用量d(噸)由下表給出。目前有兩個臨時料場位于A(5,1)，B(2,7)，日儲量各有20噸。假設(shè)從料場到工地之間均有直線道路相連。（1）試制定每天的供應(yīng)計劃，即從A，B兩料場分別向各工地運(yùn)送多少噸水泥，使總的噸千米數(shù)最小。（2）為了進(jìn)一步減少噸千米數(shù)，打算舍棄兩個臨時料場，改建兩個新的，日儲量各為20噸，問應(yīng)建在何處，節(jié)省的噸千米數(shù)有多大？（一）、建立模型

記工地的位置為(ai，bi)，水泥日用量為di，i=1,…,6;料場位置為(xj，yj)，日儲量為ej，j=1,2；從料場j向工地i的運(yùn)送量為Xij。當(dāng)用臨時料場時決策變量為：Xij，當(dāng)不用臨時料場時決策變量為：Xij，xj，yj。（二）使用臨時料場的情形

使用兩個臨時料場A(5,1)，B(2,7).求從料場j向工地i的運(yùn)送量為Xij，在各工地用量必須滿足和各料場運(yùn)送量不超過日儲量的條件下，使總的噸千米數(shù)最小，這是線性規(guī)劃問題.線性規(guī)劃模型為：設(shè)X11=X1,X21=X2,,X31=X3,X41=X4,X51=X5,,X61=X6X12=X7,X22=X8,,X32=X9,X42=X10,X52=X11,,X62=X12

編寫程序gying1.mgying1.mcleara=[1.258.750.55.7537.25];b=[1.250.754.7556.57.75];d=[3547611];x=[52];y=[17];e=[2020];fori=1:6forj=1:2

aa(i,j)=sqrt((x(j)-a(i))^2+(y(j)-b(i))^2);endendCC=[aa(:,1);aa(:,2)]';A=[111111000000;000000111111];B=[20;20];Aeq=[100000100000%從第一＼二料場運(yùn)到工地一的料

010000010000001000001000000100000100

000010000010000001000001];beq=[d(1);d(2);d(3);d(4);d(5);d(6)];VLB=[000000000000];VUB=[];x0=[123010010101];[xx,fval]=linprog(CC,A,B,Aeq,beq,VLB,VUB,x0)計算結(jié)果為：x=[3.00005.00000.00007.00000.00001.00000.00000.00004.00000.00006.000010.0000]’fval=136.2275（三）改建兩個新料場的情形

改建兩個新料場，要同時確定料場的位置(xj,yj)和運(yùn)送量Xij，在同樣條件下使總噸千米數(shù)最小。這是非線性規(guī)劃問題。非線性規(guī)劃模型為：

用lingo解此題MODEL:設(shè)X11=X1,X21=X2,,X31=X3,X41=X4,X51=X5,,X61=X6X12=X7,X22=X8,,X32=X9,X42=X10,X52=X11,,X62=X12

x1=X13,y1=X14,x2=X15,y2=X16

（1）先編寫M文件liaoch.m定義目標(biāo)函數(shù)。(2)取初值為線性規(guī)劃的計算結(jié)果及臨時料場的坐標(biāo):x0=[35070100406105127]';編寫主程序gying2.m.functionf=liaoch(x)

a=[1.258.750.55.7537.25];

b=[1.250.754.7556.57.75];

d=[3547611];

e=[2020];

f1=0;

fori=1:6

s(i)=sqrt((x(13)-a(i))^2+(x(14)-b(i))^2);

f1=s(i)*x(i)+f1;

end

f2=0;

fori=7:12

s(i)=sqrt((x(15)-a(i-6))^2+(x(16)-b(i-6))^2);

f2=s(i)*x(i)+f2;

end

f=f1+f2;

clear%x0=[35070100406105127]';%x0=[3.00005.00000.07077.000000.9293003.929306.000010.07076.38754.39435.75117.1867]';%x0=[3.00005.00000.30947.00000.01080.6798003.690605.989210.32025.53694.91945.82917.2852]';x0=[35471000005115.63484.86877.24797.7499]';A=[11111100000000000000001111110000];B=[20;20];原程序gying2.mAeq=[100000100000000001000001000000000010000010000000000100000100000000001000001000000000010000010000];beq=[3547611]';vlb=[zeros(12,1);-inf;-inf;-inf;-inf];vub=[];[x,fval,exitflag]=fmincon('liaoch',x0,A,B,Aeq,beq,vlb,vub)(3)計算結(jié)果為：x=[3.00005.00000.07077.000000.9293003.929306.000010.07076.38754.39435.75117.1867]’fval=105.4626exitflag=1(4)若修改主程序gying2.m,取初值為上面的計算結(jié)果:x0=[3.00005.00000.07077.000000.9293003.929306.000010.07076.38754.39435.75117.1867]’得結(jié)果為:x=[3.00005.00000.30947.00000.01080.6798003.690605.989210.32025.53694.91945.82917.2852]’fval=103.4760exitflag=1總的噸千米數(shù)比上面結(jié)果略優(yōu).(5)若再取剛得出的結(jié)果為初值,卻計算不出最優(yōu)解.(6)若取初值為:x0=[35471000005115.63484.86877.24797.7499]',

則計算結(jié)果為:x=[3.00005.00004.00007.00001.0000000005.000011.00005.69594.92857.25007.7500]’fval=89.8835exitflag=1總的噸千米數(shù)89.8835比上面結(jié)果更好.

gying2.m通過此例可看出fmincon函數(shù)在選取初值上的重要性.用lingo解此題的方法:Lingo的使用－－數(shù)學(xué)函數(shù)LINGO提供了大量的標(biāo)準(zhǔn)數(shù)學(xué)函數(shù)：@abs(x)返回x的絕對值@sin(x)返回x的正弦值，x采用弧度制@cos(x)返回x的余弦值@tan(x)返回x的正切值@exp(x)返回常數(shù)e的x次方@log(x)返回x的自然對數(shù)@lgm(x)返回x的gamma函數(shù)的自然對數(shù)@sign(x)如果x<0返回-1；否則，返回1@floor(x)返回x的整數(shù)部分。當(dāng)x>=0時，返回不超過x的最大整數(shù)；當(dāng)x<0時，返回不低于x的最大整數(shù)。@smax(x1,x2,…,xn)返回x1，x2，…，xn中的最大值@smin(x1,x2,…,xn)返回x1，x2，…，xn中的最小值例4.3

給定一個直角三角形，求包含該三角形的最小正方形。

LINGO代碼如下：model:Titletrangle;sets:object/1..3/:f;endsetsdata:a,b=3,4;!兩個直角邊長，修改很方便;Enddataf(1)=a*@sin(x);f(2)=b*@cos(x);f(3)=a*@cos(x)+b*@sin(x);min=@smax(f(1),f(2),f(3));@bnd(0,x,1.57);end例題:

一奶制品加工廠用牛奶生產(chǎn)A1,A2兩種奶制品，1桶牛奶可以在甲車間用12小時加工成3公斤A1，或者在乙車間用8小時加工成4公斤A2。根據(jù)市場需求，生產(chǎn)的A1,A2全部能售出，且每公斤A1獲利24元，每公斤A2獲利16元。現(xiàn)在加工廠每天能得到50桶牛奶的供應(yīng)，每天正式工人總的勞動時間480小時，并且甲車間每天至多能加工100公斤A1，乙車間的加工能力沒有限制。試為該廠制訂一個生產(chǎn)計劃，使每天獲利最大，并進(jìn)一步討論以下3個附加問題：1）若用35元可以買到1桶牛奶，應(yīng)否作這項投資？若投資，每天最多購買多少桶牛奶？2）若可以聘用臨時工人以增加勞動時間，付給臨時工人的工資最多是每小時幾元？3）由于市場需求變化，每公斤A1的獲利增加到30元，應(yīng)否改變生產(chǎn)計劃？模型代碼如下：max=72*x1+64*x2;x1+x2<=50;12*x1+8*x2<=480;3*x1<=100;求解這個模型,結(jié)果如下:Globaloptimalsolutionfoundatiteration:0Objectivevalue:3360.000

VariableValueReducedCostX120.000000.000000X230.000000.000000RowSlackorSurplusDualPrice13360.0001.00000020.00000048.0000030.0000002.000000440.000000.000000Rangesinwhichthebasisisunchanged:

ObjectiveCoefficientRanges

CurrentAllowableAllowable

VariableCoefficientIncreaseDecrease

X172.0000024.000008.000000

X264.000008.00000016.00000

RighthandSideRanges

RowCurrentAllowableAllowable

RHSIncreaseDecrease

250.0000010.000006.666667

3480.000053.3333380.00000

4100.0000INFINITY40.00000做靈敏性分析，結(jié)果如下:用lingo解此題

MODEL:

TitleLocationProblem;

sets:

demand/1..6/:a,b,d;

supply/1..2/:x,y,e;

link(demand,supply):c;

endsets

data:

!locationsforthedemand(需求點的位置);

a=1.25,8.75,0.5,5.75,3,7.25;

b=1.25,0.75,4.75,5,6.5,7.75;

!quantitiesofthedemandandsupply（供需量）;

d=3,5,4,7,6,11;e=20,20;

enddata

init:

!initiallocationsforthesupply（初始點）;

x,y=5,1,2,7;

endinit

!Objectivefunction（目標(biāo)）;

[OBJ]min=@sum(link(i,j):c(i,j)*((x(j)-a(i))^2+(y(j)-b(i))^2)^(1/2));

!demandconstraints（需求約束）;

@for(demand(i):[DEMAND_CON]@sum(supply(j):c(i,j))=d(i););

!supplyconstraints（供應(yīng)約束）;

@for(supply(i):[SUPPLY_CON]@sum(demand(j):c(j,i))<=e(i););

@for(supply:@bnd(0.5,X,8.75);@bnd(0.75,Y,7.75););

END

運(yùn)行的結(jié)果為:X(1)5.695967X(2)7.249997Y(1)4.928560Y(2)7.750000Feasiblesolutionfound.Objectivevalue:89.88349Totalsolveriterations:81

1.3應(yīng)用舉例(二)---整數(shù)規(guī)化解先看有多少種裁料方案，再進(jìn)行組合和選擇。方案：

例1合理利用線材問題現(xiàn)要做一百套鋼管,每套要長為2.9m、2.1m和1.5m的鋼管各一根。已知原料長7.4m，問應(yīng)如何下料，使用的原料最省。

設(shè)用方案Ⅰ,Ⅱ,…,Ⅷ分別裁原料鋼管x1,x2,…,x8根,則：Minz=x1+x2+x3+x4+x5+x6+x7+x82x1+x2+x3+x4≥

1002x2+x3+3x5+2x6+x7≥100x1+x3+3x4+2x6+3x7+4x8≥100x1,x2,x3,x4,x5,x6,x7,x8≥0

上題中,如果每套1.8m的鋼管要70根,要求使用的切割模式不超過3種.問應(yīng)如何下料，使用的原料最省。

解:設(shè)xi—按第i種模式切割的原料鋼管根數(shù)(i=1,2,3)r1i,r2i,r3i,r4i—第i種切割模式下,每根原料生產(chǎn)2.9m、2.1m和1.5m,1.8m的鋼管的數(shù)量.目標(biāo)函數(shù)(總根數(shù))

(套數(shù)約束)按第i種切割模式下,每根鋼管的長度限制:這是一個整數(shù)非線性規(guī)劃模型，用LINGO要運(yùn)行很長時間也難以得到最優(yōu)解。因三種切割模式的排列順序無關(guān)緊要，所以不妨增加以下約束：又鋼管的總根數(shù)有明顯的上界和下界。首先，原料鋼管的總根數(shù)不可能少于其次，考慮第一種切割模式下只生產(chǎn)2.9m鋼管，一根原料鋼管切割成2根2.9m鋼管，100套鋼管需要50根原料；如此可計算出鋼管的上界：

所以可以增加以下約束：Lingo源程序model:title鋼管下料；min=x1+x2+x3;x1*r11+x2*r12+x3*r13>=100;x1*r21+x2*r22+x3*r23>=100;x1*r31+x2*r32+x3*r33>=100;x1*r41+x2*r42+x3*r43>=70;x1+x2+x3>=105;x1+x2+x3<=152;2.9*r11+2.1*r21+1.5*r31+1.8*r41>=5.9;2.9*r12+2.1*r22+1.5*r32+1.8*r42>=5.9;2.9*r13+2.1*r23+1.5*r33+1.8*r43>=5.9;2.9*r11+2.1*r21+1.5*r31+1.8*r41<=7.4;2.9*r12+2.1*r22+1.5*r32+1.8*r42<=7.4;2.9*r13+2.1*r23+1.5*r33+1.8*r43<=7.4;x1>=x2;x2>=x3;@gin(x1);@gin(x2);@gin(x3);@gin(r11);@gin(r21);@gin(r31);@gin(r41);@gin(r12);@gin(r22);@gin(r32);@gin(r42);@gin(r13);@gin(r23);@gin(r33);@gin(r43);end運(yùn)行結(jié)果如下：Feasiblesolutionfoundatiteration:26565

ModelTitle:鋼管下料；

VariableValue

X172.00000

R111.000000

X229.00000

R120.000000

X314.00000

R132.000000

R211.000000

R221.000000

R230.000000

R310.000000

R323.000000

R331.000000

R411.000000

R420.000000

R430.000000Lingo源程序二model:Title鋼管下料-最小化鋼管根數(shù)的LINGO模型;SETS:NEEDS/1..4/:LENGTH,NUM; !定義基本集合NEEDS及其屬性LENGTH,NUM;CUTS/1..3/:X; !定義基本集合CUTS及其屬性X;PATTERNS(NEEDS,CUTS):R;!定義派生集合PATTERNS（這是一個稠密集合）及其屬性R;ENDSETSDATA: LENGTH=2.92.11.51.8; NUM=100100

10070; CAPACITY=7.4;ENDDATAmin=@SUM(CUTS(I):X(I)); !目標(biāo)函數(shù);Lingo源程序二@FOR(NEEDS(I):@SUM(CUTS(J):X(J)*R(I,J))>NUM(I));!滿足需求約束;@FOR(CUTS(J):@SUM(NEEDS(I):LENGTH(I)*R(I,J))<CAPACITY);!合理切割模式約束;@FOR(CUTS(J):@SUM(NEEDS(I):LENGTH(I)*R(I,J))>CAPACITY-@MIN(NEEDS(I):LENGTH(I))); !合理切割模式約束;@SUM(CUTS(I):X(I))>105;@SUM(CUTS(I):X(I))<115;!人為增加約束;@FOR(CUTS(I)|I#LT#@SIZE(CUTS):X(I)>X(I+1)); !人為增加約束;@FOR(CUTS(J):@GIN(X(J)));@FOR(PATTERNS(I,J):@GIN(R(I,J)));end運(yùn)行結(jié)果如下：Localoptimalsolutionfoundatiteration:64408

Objectivevalue:113.0000

ModelTitle:鋼管下料-最小化鋼管根數(shù)的LINGO模型

VariableValueReducedCost

CAPACITY