数学建模与数学实验非线性规划实验目的实验内容
2,掌握用数学软件求解优化问题．
1,直观了解非线性规划的基本内容．
1．非线性规划的基本理论．
4．实验作业．
2,用数学软件求解非线性规划．
3,钢管订购及运输优化模型．
*非线性规划的基本解法非线性规划的基本概念非线性规划返回定义 如果目标函数或约束条件中至少有一个是非线性函数,
则最优化问题就叫做 非线性规划问题,
非现性规划的基本概念一般形式,
（ 1）
其中,是定义在 Rn 上的实值函数，简记,
Xfm in
ji hgf,,
其它情况,求目标函数的最大值，或约束条件小于等于零两种情况，都可通过取其相反数化为上述一般形式．
1n
j
1n
i
1n R:h,R:g,R,RRRf
nTn RxxxX?=,,,21 L
==
=?
.,...,2,10
m;1,2,...,0..
ljXh
iXgts
j
i
定义 1 把满足问题（ 1）中条件的解 称为 可行解 （或可行点 ）,所有可行点的集合称为 可行集 （或 可行域 ）,记为 D．即问题 (1)可简记为,
XfDX?m in
定义 2 对于问题 (1),设,若存在,使得对一切
,且,都有,则称 X*是 f(X)在 D上的局部极小值点 （ 局部最优解 ）,特别地,当 时，若
,则称 X*是 f(X)在 D上的 严格局部极小值点 （ 严格局部最优解 ）,
DX?* 0
DX?
*XX
*XX?
XfXf?*
XfXf?*
定义 3 对于问题 (1),设,若对任意的,都有则称 X*是 f(X)在 D上的 全局极小值点 （ 全局最优解 ）,特别地,当时，若,则称 X*是 f(X)在 D上的 严格全局极小值点 （ 严格全局最优解 ）,
DX?* DX?
*XXXfXf?*
返回
)( nRX?
{ }nji RXXhXgXD?=?=,0,0|
,XfXf?*
非线性规划的基本解法
SUTM外点法
SUTM内点法（障碍罚函数法）
1,罚函数法
2,近似规划法返回罚函数法罚函数法 基本思想是通过构造罚函数把约束问题转化为一系列无约束最优化问题，
进而用无约束最优化方法去求解．这类方法称为 序列无约束最小化方法,简称为 SUMT
法,
其一为 SUMT外点法,其二为 SUMT内点法,
)2(,0m i n,
1
2
1
2
==
=
l
j
j
m
i
i XhMXgMXfMXT可设：
R1 m in,( 3 )nX T X M?将 问 题 （ ） 转 化 为 无 约 束 问 题,
其中 T(X,M)称为 罚函数,M称为 罚因子,带 M的项称为 罚项,
这里的罚函数只对不满足约束条件的点实行惩罚,当 时，满足各,故罚项为 0，不受惩罚．当 时,必有约束条件,故罚项大于 0，要受惩罚．
DX?
0,0 =? XhXg ii DX?
00 XhXg ii 或
SUTM外点法
m i n
0 1,2,.,,,;
s,t,( 1 )
0 1,2,.,,,.
i
j
fX
g X i m
h X j l
=
==
对 一 般 的 非 线 性 规 划,
罚函数法的 缺点,每个近似最优解 Xk往往不是容许解，而只能近似满足约束，在实际问题中这种结果可能不能使用；在解一系列无约束问题中，计算量太大，特别是随着 Mk的增大，
可能导致错误．
1,任意给定初始点 X0，取 M1>1，给定允许误差,令 k=1；
2,求无约束极值问题 的最优解，设 Xk=X(Mk)，即；
3,若存在,使,则取 Mk>M( ),
令 k=k+1返回（ 2），否则，停止迭代．得最优解,
计算时也可将收敛性判别准则 改为,
0
Rm i n,nX T X M?
Rm in,(,)n k kX T X M T X M? =
mii1 ki Xg 10,1 == MM k
0,0m i n
1
2
=
m
i
i XgM
kXX?*
ki Xg
SUTM外点法 (罚函数法 )的 迭代步骤
m i n
( 1 )
s,t,0 1,2,.,,,i
fX
g X i m?=
考 虑 问 题,
{ }00| 0,1,2,,iD X g X i m D?=? =?设 集 合,是可 行 域 中 所 有 严 格 内 点 的 集 合,
0
1
m i n,
k
kk
XD
I X r X r
这 样 问 题 （ ） 就 转 化 为 求 一 系 列 极 值 问 题,
得 （ ）,
SUTM内点法（ 障碍函数法 ）
为障碍因子,为障碍项，或其中称或：
构造障碍函数
rXgrXgr
XgrXfrXIXgrXfrXIrXI
m
i i
m
i
i
m
i i
m
i
i
==
==
=?=
11
11
1ln
1)(),(ln,,
内点法的迭代步骤
( 1 ) 给定允许误差 0,取 10,01r ；
(2) 求出约束集合 D 的一个内点 00 DX?,令 1=k ； (3) 以 01 DX k 为初始点，求解
kDX rXI,m i n 0?
，其中 0DX? 的最优解 设为 0DrXX kk?= ；
(4) 检验是否满足
=
m
i
k
i Xgr
1
ln 或
=
m
i
i
k
Xg
r
1
1
，若满足，停止迭代,令 kXX?* ；否则取
kk rr 1?=?
，令 1?= kk,
返回 (3 ）,
近似规划法的基本思想,将问题 (3)中的目标函数和约束条件近似为线性函数，并对变量的取值范围加以限制，从而得到一个近似线性规划问题，再用单纯形法求解之，
把其符合原始条件的最优解作为 (3)的解的近似．
Xf
0 ( 1,.,,,) ; 0 ( 1,,)ijg X i m h X j l? = = =
近似规划法每得到一个近似解，都从这点出发，重复以上步骤．
这样，通过求解一系列线性规划问题，产生一个由线性规划最优解组成的序列，经验表明，这样的序列往往收敛于非线性规划问题的解．
近似规划法的 算法步骤如下：
(2) 在点 kX 处，将Xf,
XhXg ji,
按泰勒级数展开并取一阶近似，得到近似线性规划问题,
T
m i n
k k k
f X f X f X X X
0 1,,
T
k k k
i i i
g X g X g X X X i m =
T
0 j 1,,
k k k
j j j
h X h X h X X X l = =；
(1) 给定初始可行点 { }112111,,,nxxxX L=,步长限制njj,,11 L=?,
步长缩小系数1,0,允许误差? > 0,令 k =1 ；
5) 判断精度：若njkj,,1 L=,则点 1?kX 为近似最优解；
否则，令njkjkj,,1 1 L==,k = k +1,返回步骤 ( 2),
( 3 ) 在上述近似线性规划问题的基础上增加一组限制步长的线性约束条件．因为线性近似通常只在展开点附近近似程度较高，故需要对变量的取值范围加以限制，所增加的约束条件是：
njxx
k
j
k
jj
,,1 L=
求解该线性规划问题，得到最优解
1?k
X ；
(4) 检验 1?kX 对原约束是否可行,若 1?kX 对原约束可行，则转步骤 (5) ；否则，缩小步长限制，令njk
j
k
j,,1 L==
，返回步骤 (3),重解当前的线性规划问题；
返回用 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．二次规划标准型为,
m i n Z =
2
1
X
T
HX + c
T
X
s.t,AX ≤ b
b e qXA e q =?
V L B ≤ X ≤ V U B
例 1 min f(x1,x2)=-2x1-6x2+x12-2x1x2+2x22
s.t,x1+x2≤2
-x1+2x2≤2
x1≥0,x2≥0
MATLAB（ youh1）
1,写成标准形式,
2,输入命令,
H=[1 -1; -1 2];
c=[-2 ;-6];A=[1 1; -1 2];b=[2;2];
Aeq=[];beq=[]; VLB=[0;0];VUB=[];
[x,z]=quadprog(H,c,A,b,Aeq,beq,VLB,VUB)
3,运算结果 为：
x =0,6667 1,3333 z = -8,2222
T
11
12
22
1 - 1 2m in (,)
1 2 6
xxz x x=?
1
2
1
2
1 1 2
1 2 2
0
0
x
x
x
x
s.t.
1,首先建立 M文件 fun.m,用来定义目标函数 F（ X）,
function f=fun(X);
f=F(X);
2．一般非线性规划标准型为,
m i n F ( X )
s,t,AX? b b e qXA e q =? G ( X ) 0?
Ceq ( X )= 0 VL B? X? V U B
其中 X为 n维变元向量,G(X)与 Ceq(X)均为非线性函数组成的向量，其他变量的含义与线性规划、二次规划中相同．用
MATLAB求解上述问题，基本步骤分三步：
2,若约束条件中有非线性约束,G (X) 0? 或 Ceq(X)= 0,
则建立 M 文件 nonlcon,m 定义函数 G(X) 与 Ceq (X),
function [G,Ceq]=nonlcon(X)
G=?
Ceq=?
3,建立主程序,求解非线性规划的函数是 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文件 迭代的初值 参数说明变量上下限注意：
[1] fmincon函数提供了大型优化算法和中型优化算法．默认时,若在 fun函数中提供了梯度（ options参数的 GradObj设置为’ on’），并且只有上下界存在或只有等式约束,fmincon函数将选择大型算法．当既有等式约束又有梯度约束时，使用中型算法．
[2] fmincon函数的中型算法使用的是序列二次规划法．在每一步迭代中求解二次规划子问题，并用 BFGS法更新拉格朗日
Hesse矩阵．
[3] fmincon函数可能会给出局部最优解，这与初值 X0的选取有关．
1,写成标准形式,
s.t.
0
0
54
632
21
21
xx
xx
2
1
0
0
x
x
2
2
2
121 2
1
2
12m i n xxxxf=
2
2
2
121 2
1
2
12m i n xxxxf=
2x1+3x2 6
s.t,x1+4x2 5
x1,x2 0
例 2
2,先建立 M-文件 fun3,m:
function f=fun3(x);
f=-x(1)-2*x(2)+(1/2)*x(1)^2+(1/2)*x(2)^2
MATLAB(youh2)
3．再建立主程序 youh2,m：
x0=[1;1];
A=[2 3 ;1 4]; b=[6;5];
Aeq=[];beq=[];
VLB=[0;0]; VUB=[];
[x,fval]=fmincon('fun3',x0,A,b,Aeq,beq,VLB,VUB)
4,运算结果为：
x = 0,7647 1,0588
fval = -2,0294
1,先建立 M文件 fun4,m定义目标函数,
function f=fun4(x);
f=exp(x(1))
*(4*x(1)^2+2*x(2)^2+4*x(1)*x(2)+2*x(2)+1);
1 221 2 1 2 2( ) e ( 4 2 4 2 1 )xf x x x x x x=
x1+x2=0
s.t,1.5+x1x2 - x1 - x2 0
-x1x2 –10 0?
例 3
2．再建立 M文件 mycon,m定义非线性约束：
function [g,ceq]=mycon(x)
g=[x(1)+x(2);1,5+x(1)*x(2)-x(1)-x(2);
-x(1)*x(2)-10];
3．主程序 youh3,m为,
x0=[-1;1];
A=[];b=[];
Aeq=[1 1];beq=[0];
vlb=[];vub=[];
[x,fval]=fmincon('fun4',x0,A,b,Aeq,beq,vlb,
vub,'mycon')
MATLAB(youh3)
4,运算结果为,
x = -1,2250 1,2250
fval = 1,8951
例 4
12
22
1 1 2
22
2 1 2
12
m in 2
s,t,2 5 0
7 0
0 5,0 10
f X x x
g X x x
g X x x
xx
=
=
=
1．先建立 M文件 fun,m定义目标函数,
function f=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=[0 0];VUB=[5 10];
[x,fval,exitflag,output]
=fmincon('fun',x0,[],[],[],[],
VLB,VUB,'mycon2')
MATLAB(fxx(fun))
4,运算结果为,
x =
4,0000
3,0000
fval =-11,0000
exitflag = 1
output =
iterations,4
funcCount,17
stepsize,1
algorithm,[1x44 char]
firstorderopt,[]
cgiterations,[]
返回应用实例,供应与选址某公司有 6个建筑工地要开工，每个工地的位置（用平面坐标系
a,b表示，距离单位,km）及水泥日用量 d(t)由下表给出．目前有两个临时料场位于 A(5,1),B(2,7)，日储量各有 20t．假设从料场到工地之间均有直线道路相连．
（ 1）试制定每天的供应计划，即从 A,B两料场分别向各工地运送多少水泥，可使总的吨千米数最小．
（ 2）为了进一步减少吨千米数，打算舍弃两个临时料场，改建两个新的，日储量各为 20t，问应建在何处，节省的吨千米数有多大？
工地位置（ a,b ）及水泥日用量 d
1 2 3 4 5 6
a 1.25 8.75 0.5 5.75 3 7.25
b 1.25 0.75 4.75 5 6.5 7.25
d 3 5 4 7 6 11
（一）建立模型记工地的位置为 (ai,bi)，水泥日用量为 di,i=1,…,6; 料场位置为
(xj,yj)，日储量为 ej,j=1,2；料场 j向工地 i的运送量为 Xij．目标函数为,
= =
=
2
1
6
1
22 )()(m i n
j i
ijijij byaxXf
约束条件为,
2,1,
6,,2,1,
6
1
2
1
=?
==
=
=
jeX
idX
j
i
ij
i
j
ij
L
当用临时料场时决策变量为,Xij，
当不用临时料场时决策变量为,Xij,xj,yj．
（二）使用临时料场的情形使用两个临时料场 A(5,1),B(2,7)．求从料场 j向工地 i的运送量 Xij,在各工地用量必须满足和各料场运送量不超过日储量的条件下，使总的吨千米数最小，这是线性规划问题,线性规划模型为：
= =
=
2
1
6
1
),(m in
j i
ij
Xjiaaf
2,1,
6,,2,1,s,t,
6
1
2
1
=?
==
=
=
jeX
idX
j
i
ij
i
j
ij
L
其中 22 )()(),( ijij byaxjiaa=,i = 1,2,…,6,j = 1,2,为常数,
设 X11=X1,X21= X 2,,X31= X 3,X41= X 4,X51= X 5,,X61= X 6
X12= X 7,X22= X 8,,X32= X 9,X42= X 10,X52= X 11,,X62= X 12
编写程序 gying1,m MATLAB（ gying1）
计算结果为：
x =[ 3,0000 5,0000 0,0000 7,0000 0,0000 1,0000 0,0000
0,0000 4,0000 0,0000 6,0000 10,0000]’
fval = 136,2275
即由料场 A,B 向 6 个工地运料方案为,
1 2 3 4 5 6
料场 A 3 5 0 7 0 1
料场 B 0 0 4 0 6 10
总的吨千米数为 13 6.2 275,
（三）改建两个新料场的情形改建两个新料场，要同时确定料场的位置 (xj,yj)和运送量
Xij，在同样条件下使总吨千米数最小．这是非线性规划问题．非线性规划模型为：
= =
=
2
1
6
1
22 )()(m i n
j i
ijijij byaxXf
2
1
6
1
s,t.,1,2,,6
,1,2
ij i
j
ij j
i
X d i
X e j
=
=
==
=
设 X11=X1,X21= X 2,X31= X 3,X41= X 4,X51= X 5,,X61= X 6
X12= X 7,X22= X 8,X32= X 9,X42= X 10,X52= X 11,X62= X 12
x1=X13,y1=X14,x2=X15,y2=X16
（ 1）先编写 M文件 liaoch,m定义目标函数．
MATLAB（ liaoch）
(2) 取初值为线性规划的计算结果及临时料场的坐标,
x0=[3 5 0 7 0 1 0 0 4 0 6 10 5 1 2 7]';
编写主程序 gying2,m．
MATLAB（ gying2）
(3) 计算结果为：
x=[ 3,0000 5,0000 0,0707 7,0000 0 0,9293
0 0 3,9293 0 6,0000 10,0707
6,3875 4,3943 5,7511 7,1867]’
fval = 105,4626
exitflag = 1 即两个新料场的坐标分别为 ( 6.38 7 5,4,394 3),( 5,75 1 1,7.18 67),由料场 A,B 向 6 个工地运料方案为,
1 2 3 4 5 6
料场 A 3 5 0.0707 7 0 0.9293
料场 B 0 0 3.9293 0 6 10.0 70 7
总的吨千米数为 105,462 6,比用临时料场节省的 吨千米 数 31,
(4) 若修改主程序 gying2,m,取初值为上面的计算结果,
x0=[ 3,0000 5,0000 0,0707 7,0000 0
0,9293 0 0 3,9293 0 6,0000 10,0707
6,3875 4,3943 5,7511 7,1867]’
则得结果为,
x=[3,0000 5,0000 0,3094 7,0000
0,0108 0,6798 0 0 3,6906 0
5,9892 10,3202 5,5369 4,9194
5,8291 7,2852]’
fval =103,4760
exitflag = 1
总的吨千米数比上面结果略优．
(5) 若再取刚得出的结果为初值,却计算不出最优解．
MATLAB（ gying2）
MATLAB（ gying2）
(6) 若取初值为,
x0=[3 5 4 7 1 0 0 0 0 0 5 11 5,6348 4,8687 7,2479 7,7499]',
则计算结果为,
x=[3,0000 5,0000 4,0000 7,0000 1,0000 0 0 0 0 0
5,0000 11,0000 5,6959 4,9285 7,2500 7,7500]’
fval =89,8835
exitflag = 1
总的吨千米数 89,8835比上面结果更好．
通过此例可看出 fmincon函数在选取初值上的重要性．
MATLAB（ gying2）
返回钢管订购及运输优化模型
2000年“网易杯”全国大学生数学建模竞赛 B题符号说明：
的距离到：表示结点 11, jjjj AAA
个结点的钢管数量个工厂运到第：表示从第 jix ij
所需要的钢管数量：表示结点 jj An
个结点的最大生产能力：表示第 js j
jA 1?jA
个结点的运费和成本个工厂运到第：单位产品从第 jia ij
是已知的量，是待求的变量，而其中 1.,,,?jjijijjij Asnatx
)0( 1 =t之间的平衡点到：表示结点 1?jjj AAt
jt
20
))(1()1(
)
2
))(1(
2
)1(
(1.0
))}(21()21{(1.0
1.1.
1.1.
1.1.
jjjjjjjj
jjjjjjjj
jjjjjj
tAtAtt
tAtAtt
tAtC
=
=
=
LL
1．铺设总费用：
=
=
14
1
1.
j
jjCC
2．成本及运输总费用：
= =
=
7
1
15
2i j
ijij axW
总费用 =铺设总费用 +成本及运输总费用 =C+W
模型的分析与建立
1 4 1 5 7
.1
1 2 1
7
1
1 5 1 5
22
.1
j 2,3,,1 5
5 0 0 0
s.t,
0 1,,7,2,,1 5
0
j j ij ij
j j i
ij j
i
ij i ij
jj
ij
j j j
f C x a
xn
x s x
x i j
tA
= = =
=
==
=
==
=
= =
或建立模型模型求解利用 MATLAB软件包求解得：
钢 厂 S 1 S 2 S 3 S 4 S 5 S 6 S 7
订购量 800 800 1000 0 1015 1550 0
总费用 12786 32
订购量 A
2
A
3
A
4
A
5
A
6
A
7
A
8
A
9
A
10
A
11
A
12
A
13
A
14
A
15
S
1
800 0 201 133 200 266 0 0 0 0 0 0 0 0 0
S
2
800 179 11 14 295 0 0 300 0 0 0 0 0 0 0
S
3
1000 139 11 186 0 0 0 664 0 0 0 0 0 0 0
S
4
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
S
5
1015 0 358 242 0 0 0 0 0 0 415 0 0 0 0
S
6
1556 0 0 0 0 0 0 0 0 0 351 86 333 621 165
S
7
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
订购和运输方案表返回某厂向用户提供发动机，合同规定，第一、二、
三季度末分别交货 40台,60台,80台．每季度的生产费用为 （单位,元）,其中 x是该季度生产的台数．若交货后有剩余，可用于下季度交货，
但需支付存储费，每台每季度 c元．已知工厂每季度最大生产能力为 100台，第一季度开始时无存货，设
a=50,b=0.2,c=4，问,工厂应如何安排生产计划，
才能既满足合同又使总费用最低．讨论 a,b,c变化对计划的影响，并作出合理的解释．
练习 1
2bxaxxf?=
一基金管理人的工作是,每天将现有的美元、英镑、马克和 日元四种货币按当天汇率相互兑换，使在满足需要的条件下，按美元计算的价值最高．设某天的汇率、现有货币和当天需求如下,
美元 英镑 马克 日元现有量
)10(
8
需求量
)10(
8
美元 1 0,58928 1,743 138,3 8 6
英镑 1,697 1 2,9579 234,7 1 3
马克 0,57372 0,33808 1 79,346 8 1
日元 0,00 7233 0,00426 0,0126 1 0 10
问该天基金管理人应如何操作,(,按美元计算的价值,
指 兑 入,兑 出 汇 率 的 平 均 值,如 1 英 镑 相 当 于
2
58928.01697.1?
=1.6969 93 美元,)
练习 2
返回谢 谢！
2,掌握用数学软件求解优化问题．
1,直观了解非线性规划的基本内容．
1．非线性规划的基本理论．
4．实验作业．
2,用数学软件求解非线性规划．
3,钢管订购及运输优化模型．
*非线性规划的基本解法非线性规划的基本概念非线性规划返回定义 如果目标函数或约束条件中至少有一个是非线性函数,
则最优化问题就叫做 非线性规划问题,
非现性规划的基本概念一般形式,
（ 1）
其中,是定义在 Rn 上的实值函数，简记,
Xfm in
ji hgf,,
其它情况,求目标函数的最大值，或约束条件小于等于零两种情况，都可通过取其相反数化为上述一般形式．
1n
j
1n
i
1n R:h,R:g,R,RRRf
nTn RxxxX?=,,,21 L
==
=?
.,...,2,10
m;1,2,...,0..
ljXh
iXgts
j
i
定义 1 把满足问题（ 1）中条件的解 称为 可行解 （或可行点 ）,所有可行点的集合称为 可行集 （或 可行域 ）,记为 D．即问题 (1)可简记为,
XfDX?m in
定义 2 对于问题 (1),设,若存在,使得对一切
,且,都有,则称 X*是 f(X)在 D上的局部极小值点 （ 局部最优解 ）,特别地,当 时，若
,则称 X*是 f(X)在 D上的 严格局部极小值点 （ 严格局部最优解 ）,
DX?* 0
DX?
*XX
*XX?
XfXf?*
XfXf?*
定义 3 对于问题 (1),设,若对任意的,都有则称 X*是 f(X)在 D上的 全局极小值点 （ 全局最优解 ）,特别地,当时，若,则称 X*是 f(X)在 D上的 严格全局极小值点 （ 严格全局最优解 ）,
DX?* DX?
*XXXfXf?*
返回
)( nRX?
{ }nji RXXhXgXD?=?=,0,0|
,XfXf?*
非线性规划的基本解法
SUTM外点法
SUTM内点法（障碍罚函数法）
1,罚函数法
2,近似规划法返回罚函数法罚函数法 基本思想是通过构造罚函数把约束问题转化为一系列无约束最优化问题，
进而用无约束最优化方法去求解．这类方法称为 序列无约束最小化方法,简称为 SUMT
法,
其一为 SUMT外点法,其二为 SUMT内点法,
)2(,0m i n,
1
2
1
2
==
=
l
j
j
m
i
i XhMXgMXfMXT可设：
R1 m in,( 3 )nX T X M?将 问 题 （ ） 转 化 为 无 约 束 问 题,
其中 T(X,M)称为 罚函数,M称为 罚因子,带 M的项称为 罚项,
这里的罚函数只对不满足约束条件的点实行惩罚,当 时，满足各,故罚项为 0，不受惩罚．当 时,必有约束条件,故罚项大于 0，要受惩罚．
DX?
0,0 =? XhXg ii DX?
00 XhXg ii 或
SUTM外点法
m i n
0 1,2,.,,,;
s,t,( 1 )
0 1,2,.,,,.
i
j
fX
g X i m
h X j l
=
==
对 一 般 的 非 线 性 规 划,
罚函数法的 缺点,每个近似最优解 Xk往往不是容许解，而只能近似满足约束，在实际问题中这种结果可能不能使用；在解一系列无约束问题中，计算量太大，特别是随着 Mk的增大，
可能导致错误．
1,任意给定初始点 X0，取 M1>1，给定允许误差,令 k=1；
2,求无约束极值问题 的最优解，设 Xk=X(Mk)，即；
3,若存在,使,则取 Mk>M( ),
令 k=k+1返回（ 2），否则，停止迭代．得最优解,
计算时也可将收敛性判别准则 改为,
0
Rm i n,nX T X M?
Rm in,(,)n k kX T X M T X M? =
mii1 ki Xg 10,1 == MM k
0,0m i n
1
2
=
m
i
i XgM
kXX?*
ki Xg
SUTM外点法 (罚函数法 )的 迭代步骤
m i n
( 1 )
s,t,0 1,2,.,,,i
fX
g X i m?=
考 虑 问 题,
{ }00| 0,1,2,,iD X g X i m D?=? =?设 集 合,是可 行 域 中 所 有 严 格 内 点 的 集 合,
0
1
m i n,
k
kk
XD
I X r X r
这 样 问 题 （ ） 就 转 化 为 求 一 系 列 极 值 问 题,
得 （ ）,
SUTM内点法（ 障碍函数法 ）
为障碍因子,为障碍项，或其中称或：
构造障碍函数
rXgrXgr
XgrXfrXIXgrXfrXIrXI
m
i i
m
i
i
m
i i
m
i
i
==
==
=?=
11
11
1ln
1)(),(ln,,
内点法的迭代步骤
( 1 ) 给定允许误差 0,取 10,01r ；
(2) 求出约束集合 D 的一个内点 00 DX?,令 1=k ； (3) 以 01 DX k 为初始点，求解
kDX rXI,m i n 0?
，其中 0DX? 的最优解 设为 0DrXX kk?= ；
(4) 检验是否满足
=
m
i
k
i Xgr
1
ln 或
=
m
i
i
k
Xg
r
1
1
，若满足，停止迭代,令 kXX?* ；否则取
kk rr 1?=?
，令 1?= kk,
返回 (3 ）,
近似规划法的基本思想,将问题 (3)中的目标函数和约束条件近似为线性函数，并对变量的取值范围加以限制，从而得到一个近似线性规划问题，再用单纯形法求解之，
把其符合原始条件的最优解作为 (3)的解的近似．
Xf
0 ( 1,.,,,) ; 0 ( 1,,)ijg X i m h X j l? = = =
近似规划法每得到一个近似解，都从这点出发，重复以上步骤．
这样，通过求解一系列线性规划问题，产生一个由线性规划最优解组成的序列，经验表明，这样的序列往往收敛于非线性规划问题的解．
近似规划法的 算法步骤如下：
(2) 在点 kX 处，将Xf,
XhXg ji,
按泰勒级数展开并取一阶近似，得到近似线性规划问题,
T
m i n
k k k
f X f X f X X X
0 1,,
T
k k k
i i i
g X g X g X X X i m =
T
0 j 1,,
k k k
j j j
h X h X h X X X l = =；
(1) 给定初始可行点 { }112111,,,nxxxX L=,步长限制njj,,11 L=?,
步长缩小系数1,0,允许误差? > 0,令 k =1 ；
5) 判断精度：若njkj,,1 L=,则点 1?kX 为近似最优解；
否则，令njkjkj,,1 1 L==,k = k +1,返回步骤 ( 2),
( 3 ) 在上述近似线性规划问题的基础上增加一组限制步长的线性约束条件．因为线性近似通常只在展开点附近近似程度较高，故需要对变量的取值范围加以限制，所增加的约束条件是：
njxx
k
j
k
jj
,,1 L=
求解该线性规划问题，得到最优解
1?k
X ；
(4) 检验 1?kX 对原约束是否可行,若 1?kX 对原约束可行，则转步骤 (5) ；否则，缩小步长限制，令njk
j
k
j,,1 L==
，返回步骤 (3),重解当前的线性规划问题；
返回用 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．二次规划标准型为,
m i n Z =
2
1
X
T
HX + c
T
X
s.t,AX ≤ b
b e qXA e q =?
V L B ≤ X ≤ V U B
例 1 min f(x1,x2)=-2x1-6x2+x12-2x1x2+2x22
s.t,x1+x2≤2
-x1+2x2≤2
x1≥0,x2≥0
MATLAB（ youh1）
1,写成标准形式,
2,输入命令,
H=[1 -1; -1 2];
c=[-2 ;-6];A=[1 1; -1 2];b=[2;2];
Aeq=[];beq=[]; VLB=[0;0];VUB=[];
[x,z]=quadprog(H,c,A,b,Aeq,beq,VLB,VUB)
3,运算结果 为：
x =0,6667 1,3333 z = -8,2222
T
11
12
22
1 - 1 2m in (,)
1 2 6
xxz x x=?
1
2
1
2
1 1 2
1 2 2
0
0
x
x
x
x
s.t.
1,首先建立 M文件 fun.m,用来定义目标函数 F（ X）,
function f=fun(X);
f=F(X);
2．一般非线性规划标准型为,
m i n F ( X )
s,t,AX? b b e qXA e q =? G ( X ) 0?
Ceq ( X )= 0 VL B? X? V U B
其中 X为 n维变元向量,G(X)与 Ceq(X)均为非线性函数组成的向量，其他变量的含义与线性规划、二次规划中相同．用
MATLAB求解上述问题，基本步骤分三步：
2,若约束条件中有非线性约束,G (X) 0? 或 Ceq(X)= 0,
则建立 M 文件 nonlcon,m 定义函数 G(X) 与 Ceq (X),
function [G,Ceq]=nonlcon(X)
G=?
Ceq=?
3,建立主程序,求解非线性规划的函数是 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文件 迭代的初值 参数说明变量上下限注意：
[1] fmincon函数提供了大型优化算法和中型优化算法．默认时,若在 fun函数中提供了梯度（ options参数的 GradObj设置为’ on’），并且只有上下界存在或只有等式约束,fmincon函数将选择大型算法．当既有等式约束又有梯度约束时，使用中型算法．
[2] fmincon函数的中型算法使用的是序列二次规划法．在每一步迭代中求解二次规划子问题，并用 BFGS法更新拉格朗日
Hesse矩阵．
[3] fmincon函数可能会给出局部最优解，这与初值 X0的选取有关．
1,写成标准形式,
s.t.
0
0
54
632
21
21
xx
xx
2
1
0
0
x
x
2
2
2
121 2
1
2
12m i n xxxxf=
2
2
2
121 2
1
2
12m i n xxxxf=
2x1+3x2 6
s.t,x1+4x2 5
x1,x2 0
例 2
2,先建立 M-文件 fun3,m:
function f=fun3(x);
f=-x(1)-2*x(2)+(1/2)*x(1)^2+(1/2)*x(2)^2
MATLAB(youh2)
3．再建立主程序 youh2,m：
x0=[1;1];
A=[2 3 ;1 4]; b=[6;5];
Aeq=[];beq=[];
VLB=[0;0]; VUB=[];
[x,fval]=fmincon('fun3',x0,A,b,Aeq,beq,VLB,VUB)
4,运算结果为：
x = 0,7647 1,0588
fval = -2,0294
1,先建立 M文件 fun4,m定义目标函数,
function f=fun4(x);
f=exp(x(1))
*(4*x(1)^2+2*x(2)^2+4*x(1)*x(2)+2*x(2)+1);
1 221 2 1 2 2( ) e ( 4 2 4 2 1 )xf x x x x x x=
x1+x2=0
s.t,1.5+x1x2 - x1 - x2 0
-x1x2 –10 0?
例 3
2．再建立 M文件 mycon,m定义非线性约束：
function [g,ceq]=mycon(x)
g=[x(1)+x(2);1,5+x(1)*x(2)-x(1)-x(2);
-x(1)*x(2)-10];
3．主程序 youh3,m为,
x0=[-1;1];
A=[];b=[];
Aeq=[1 1];beq=[0];
vlb=[];vub=[];
[x,fval]=fmincon('fun4',x0,A,b,Aeq,beq,vlb,
vub,'mycon')
MATLAB(youh3)
4,运算结果为,
x = -1,2250 1,2250
fval = 1,8951
例 4
12
22
1 1 2
22
2 1 2
12
m in 2
s,t,2 5 0
7 0
0 5,0 10
f X x x
g X x x
g X x x
xx
=
=
=
1．先建立 M文件 fun,m定义目标函数,
function f=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=[0 0];VUB=[5 10];
[x,fval,exitflag,output]
=fmincon('fun',x0,[],[],[],[],
VLB,VUB,'mycon2')
MATLAB(fxx(fun))
4,运算结果为,
x =
4,0000
3,0000
fval =-11,0000
exitflag = 1
output =
iterations,4
funcCount,17
stepsize,1
algorithm,[1x44 char]
firstorderopt,[]
cgiterations,[]
返回应用实例,供应与选址某公司有 6个建筑工地要开工，每个工地的位置（用平面坐标系
a,b表示，距离单位,km）及水泥日用量 d(t)由下表给出．目前有两个临时料场位于 A(5,1),B(2,7)，日储量各有 20t．假设从料场到工地之间均有直线道路相连．
（ 1）试制定每天的供应计划，即从 A,B两料场分别向各工地运送多少水泥，可使总的吨千米数最小．
（ 2）为了进一步减少吨千米数，打算舍弃两个临时料场，改建两个新的，日储量各为 20t，问应建在何处，节省的吨千米数有多大？
工地位置（ a,b ）及水泥日用量 d
1 2 3 4 5 6
a 1.25 8.75 0.5 5.75 3 7.25
b 1.25 0.75 4.75 5 6.5 7.25
d 3 5 4 7 6 11
（一）建立模型记工地的位置为 (ai,bi)，水泥日用量为 di,i=1,…,6; 料场位置为
(xj,yj)，日储量为 ej,j=1,2；料场 j向工地 i的运送量为 Xij．目标函数为,
= =
=
2
1
6
1
22 )()(m i n
j i
ijijij byaxXf
约束条件为,
2,1,
6,,2,1,
6
1
2
1
=?
==
=
=
jeX
idX
j
i
ij
i
j
ij
L
当用临时料场时决策变量为,Xij，
当不用临时料场时决策变量为,Xij,xj,yj．
（二）使用临时料场的情形使用两个临时料场 A(5,1),B(2,7)．求从料场 j向工地 i的运送量 Xij,在各工地用量必须满足和各料场运送量不超过日储量的条件下，使总的吨千米数最小，这是线性规划问题,线性规划模型为：
= =
=
2
1
6
1
),(m in
j i
ij
Xjiaaf
2,1,
6,,2,1,s,t,
6
1
2
1
=?
==
=
=
jeX
idX
j
i
ij
i
j
ij
L
其中 22 )()(),( ijij byaxjiaa=,i = 1,2,…,6,j = 1,2,为常数,
设 X11=X1,X21= X 2,,X31= X 3,X41= X 4,X51= X 5,,X61= X 6
X12= X 7,X22= X 8,,X32= X 9,X42= X 10,X52= X 11,,X62= X 12
编写程序 gying1,m MATLAB（ gying1）
计算结果为：
x =[ 3,0000 5,0000 0,0000 7,0000 0,0000 1,0000 0,0000
0,0000 4,0000 0,0000 6,0000 10,0000]’
fval = 136,2275
即由料场 A,B 向 6 个工地运料方案为,
1 2 3 4 5 6
料场 A 3 5 0 7 0 1
料场 B 0 0 4 0 6 10
总的吨千米数为 13 6.2 275,
（三）改建两个新料场的情形改建两个新料场，要同时确定料场的位置 (xj,yj)和运送量
Xij，在同样条件下使总吨千米数最小．这是非线性规划问题．非线性规划模型为：
= =
=
2
1
6
1
22 )()(m i n
j i
ijijij byaxXf
2
1
6
1
s,t.,1,2,,6
,1,2
ij i
j
ij j
i
X d i
X e j
=
=
==
=
设 X11=X1,X21= X 2,X31= X 3,X41= X 4,X51= X 5,,X61= X 6
X12= X 7,X22= X 8,X32= X 9,X42= X 10,X52= X 11,X62= X 12
x1=X13,y1=X14,x2=X15,y2=X16
（ 1）先编写 M文件 liaoch,m定义目标函数．
MATLAB（ liaoch）
(2) 取初值为线性规划的计算结果及临时料场的坐标,
x0=[3 5 0 7 0 1 0 0 4 0 6 10 5 1 2 7]';
编写主程序 gying2,m．
MATLAB（ gying2）
(3) 计算结果为：
x=[ 3,0000 5,0000 0,0707 7,0000 0 0,9293
0 0 3,9293 0 6,0000 10,0707
6,3875 4,3943 5,7511 7,1867]’
fval = 105,4626
exitflag = 1 即两个新料场的坐标分别为 ( 6.38 7 5,4,394 3),( 5,75 1 1,7.18 67),由料场 A,B 向 6 个工地运料方案为,
1 2 3 4 5 6
料场 A 3 5 0.0707 7 0 0.9293
料场 B 0 0 3.9293 0 6 10.0 70 7
总的吨千米数为 105,462 6,比用临时料场节省的 吨千米 数 31,
(4) 若修改主程序 gying2,m,取初值为上面的计算结果,
x0=[ 3,0000 5,0000 0,0707 7,0000 0
0,9293 0 0 3,9293 0 6,0000 10,0707
6,3875 4,3943 5,7511 7,1867]’
则得结果为,
x=[3,0000 5,0000 0,3094 7,0000
0,0108 0,6798 0 0 3,6906 0
5,9892 10,3202 5,5369 4,9194
5,8291 7,2852]’
fval =103,4760
exitflag = 1
总的吨千米数比上面结果略优．
(5) 若再取刚得出的结果为初值,却计算不出最优解．
MATLAB（ gying2）
MATLAB（ gying2）
(6) 若取初值为,
x0=[3 5 4 7 1 0 0 0 0 0 5 11 5,6348 4,8687 7,2479 7,7499]',
则计算结果为,
x=[3,0000 5,0000 4,0000 7,0000 1,0000 0 0 0 0 0
5,0000 11,0000 5,6959 4,9285 7,2500 7,7500]’
fval =89,8835
exitflag = 1
总的吨千米数 89,8835比上面结果更好．
通过此例可看出 fmincon函数在选取初值上的重要性．
MATLAB（ gying2）
返回钢管订购及运输优化模型
2000年“网易杯”全国大学生数学建模竞赛 B题符号说明：
的距离到：表示结点 11, jjjj AAA
个结点的钢管数量个工厂运到第：表示从第 jix ij
所需要的钢管数量：表示结点 jj An
个结点的最大生产能力：表示第 js j
jA 1?jA
个结点的运费和成本个工厂运到第：单位产品从第 jia ij
是已知的量，是待求的变量，而其中 1.,,,?jjijijjij Asnatx
)0( 1 =t之间的平衡点到：表示结点 1?jjj AAt
jt
20
))(1()1(
)
2
))(1(
2
)1(
(1.0
))}(21()21{(1.0
1.1.
1.1.
1.1.
jjjjjjjj
jjjjjjjj
jjjjjj
tAtAtt
tAtAtt
tAtC
=
=
=
LL
1．铺设总费用：
=
=
14
1
1.
j
jjCC
2．成本及运输总费用：
= =
=
7
1
15
2i j
ijij axW
总费用 =铺设总费用 +成本及运输总费用 =C+W
模型的分析与建立
1 4 1 5 7
.1
1 2 1
7
1
1 5 1 5
22
.1
j 2,3,,1 5
5 0 0 0
s.t,
0 1,,7,2,,1 5
0
j j ij ij
j j i
ij j
i
ij i ij
jj
ij
j j j
f C x a
xn
x s x
x i j
tA
= = =
=
==
=
==
=
= =
或建立模型模型求解利用 MATLAB软件包求解得：
钢 厂 S 1 S 2 S 3 S 4 S 5 S 6 S 7
订购量 800 800 1000 0 1015 1550 0
总费用 12786 32
订购量 A
2
A
3
A
4
A
5
A
6
A
7
A
8
A
9
A
10
A
11
A
12
A
13
A
14
A
15
S
1
800 0 201 133 200 266 0 0 0 0 0 0 0 0 0
S
2
800 179 11 14 295 0 0 300 0 0 0 0 0 0 0
S
3
1000 139 11 186 0 0 0 664 0 0 0 0 0 0 0
S
4
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
S
5
1015 0 358 242 0 0 0 0 0 0 415 0 0 0 0
S
6
1556 0 0 0 0 0 0 0 0 0 351 86 333 621 165
S
7
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
订购和运输方案表返回某厂向用户提供发动机，合同规定，第一、二、
三季度末分别交货 40台,60台,80台．每季度的生产费用为 （单位,元）,其中 x是该季度生产的台数．若交货后有剩余，可用于下季度交货，
但需支付存储费，每台每季度 c元．已知工厂每季度最大生产能力为 100台，第一季度开始时无存货，设
a=50,b=0.2,c=4，问,工厂应如何安排生产计划，
才能既满足合同又使总费用最低．讨论 a,b,c变化对计划的影响，并作出合理的解释．
练习 1
2bxaxxf?=
一基金管理人的工作是,每天将现有的美元、英镑、马克和 日元四种货币按当天汇率相互兑换，使在满足需要的条件下，按美元计算的价值最高．设某天的汇率、现有货币和当天需求如下,
美元 英镑 马克 日元现有量
)10(
8
需求量
)10(
8
美元 1 0,58928 1,743 138,3 8 6
英镑 1,697 1 2,9579 234,7 1 3
马克 0,57372 0,33808 1 79,346 8 1
日元 0,00 7233 0,00426 0,0126 1 0 10
问该天基金管理人应如何操作,(,按美元计算的价值,
指 兑 入,兑 出 汇 率 的 平 均 值,如 1 英 镑 相 当 于
2
58928.01697.1?
=1.6969 93 美元,)
练习 2
返回谢 谢！
