数模往届2013无约束

4.

无约束规划无约束规划建模——引例无约束规划模型无约束规划的示意图无约束规划的性质

无约束规划重要算法用MATLAB解无约束规划引例1：某城市要选定一个供应中心(供电、供水、公天然气、网络接入等)的位置，该中心向城市中由固定位置的m个用户提供服务。问如何确定这个中心的位置才是最合理的？解

设用户位置为(ai,bi)(i=1,2,...,m)已知，而供应中心的位置是(x1,x2)；假设运输线路不受道路影响，只考虑直线距离，则可建立求总距离最短的数学模型为：mi

1

a

bx

x

min

f

(

x1

,

x2

)

2

21

i

2

i其中：c——运价(元/吨公里)；wi——第i个用户的需求量.w

x

a

x

b

mi

1min

f

(

x1

,

x2

)

c2

2i

1

i

2

i若进一步考虑各单位的用量wi不同，而运价c按距离计算，则可建立求总总费用最小的数学模型为：“定位问题”——厂址选择等即要求均方误差最小——最小二乘问题，包括线性最小二乘和非线性最小二乘问题，在数学建模中有非常重要的应用。mnj

121

f

(

x

,,

x

)“最小二乘问题”——解方程组引例2：在某些实际问题中，往往要求决策变量x1,x2,...,xn满足(或近似满足)某些平衡条件，它表现为要求x1,x2,...,xn满足方程组f

j

(

x1

,,

xn

)

0,

j

1,2,,

p可转化为求解问题：min无约束规划模型标准形式：X

Rnmin

f

X

(

其中

f

:

R

n

R

)max

f

X

min

f

X

XRn

XRn转化：x1x2f

(

x1

x2

)O1xx2等高线03510XX1X

2f

(

X

0

)

f

(

X1

)

f

(

X

2

)无约束规划的示意图多局部极小f

0

.

298f

0f

0

.

298唯一极小(全局极小)22

1

221

1

21

2f

(xx

3x

xx

)

2x

2x

x

2

2

f

2

fx

xx1x2n

2x

xn

1x

x2

n1

n

2

f

x

x

2

f

2

f

2

f

22

2

f

2

f

2

fx1

221x2

f

(

X

)

H

(

X

)

nx无约束规划的性质梯度x

x

xf

(

X

)

T

f

(

X

)

f

(

X

)f

(

X

)

1

2

n

,

,,

列向量梯度的特点：1)函数f(X)在X处增长最快的方向，负梯度方向是函数f(X)在X处下降最快的方向；2)梯度的模是函数最快的变化率;3)梯度为零是驻点，极值点是驻点，驻点不一定是极值点。

方、对x

x

称矩阵海森(Hessian)矩阵无约束规划的性质梯度为0向量的点可能是局部最优点，需要用海森矩阵类型进一步判定；凸函数的驻点必是全局最优点；补充1)

f(X)是凸函数对任意定义域D中的元素X、Y和任意的数0≤a≤1，有：f(aX＋(1－a)Y)

≤af(X)+

(1－a)

f(Y)补充2)

对称矩阵的分类类别A正定A半正定A负定半负定不定定义特征值>0特征值≥0(有0)特征值<0特征值≤0(有0)其它判别方法主子式都>0|A|=0;主子式≥0－A正定－A半正定无约束规划的基本算法基本思路：先求驻点，在判定是否是极值点；(初始值),寻求使函数值下降一个方向,沿此方向找到一个更好的估计值,这样往复迭代,直到某一点不能找到下降方向则终止算法.迭代算法的基本框架：选定初始点X0，令k＝0；在Xk处选定下降方向Sk(若找不到则终止算法);(3)

从Xk出发沿Sk一维搜索,找到Xk+1=Xk

+λkSk

,使得f(Xk+1)<f(Xk)；(4)

从k=k+1,转(2).X

kX

k

1X

k

2SS

k基本方法——迭代算法:先给定极小点的估计值Sk

2kS

k无约束规划的基本算法一维搜索基本原则：1)最优原则：2)

可接受点原则：

Sk

)kk

kk

kS

)

min

f

(

X

:

f

(

X

0S

k

)

f

(

X

k

)kk

k

:

f

(

X

一维搜索方法：0.618法、插值法等下降方向：与梯度点乘为负值的方向f

(

X

k

)T

Sk

0

f

(

X

S

)

f

(

X

)

f

(

X

)T

S

o(

S

)

0无约束规划的基本算法下降方向选择方法——算法分类1.

最速下降法——梯度法Sk

f

(

X

k

)2.

共轭梯度法(包含多种方法)——梯度法的修正3.

牛顿法Sk

2

f

(

X

k

)1

f

(

X

k

)f

(

X

S

)

f

(

X

)

f

(

X

)T

S

ST

2

f

(

X

)S

o(

S

2

)

0牛顿法特点：二次函数1次收敛；单步迭代计算量大，且要求海森矩阵可逆。无约束规划的基本算法拟牛顿法为克服牛顿法的缺点，同时保持较快收敛速度的优点，利用第k

步和第k+1

步得到的X

k

，X

k

1

，f

(X

k

)，f

(X

k

1

)，构造一个正定矩阵H

k

1

近似代替(2

f

(X

k

))1

，将牛顿方向改为：S

k

1

H

k

1f

(

xk

1

)从而得到下降方向.拟牛顿法特点：较梯度法收敛速度快速，较牛顿法稳定且计算量小。信赖域算法——对大型问题稳定性好，有成功的应用无约束规划的基本算法两个著名的拟牛顿迭代公式:kH

k

1

H

k

X

k

(gk

)T

H

k

H

k

gk

(X

k

)T(gk

)T

X

k(gk

)T

xk

(gk

)T

X(gk

)T

H

k

gk

xk

(X

k

)T

1

H

H

(gk

)T

X

k

(gk

)T

H

k

gkX

k

(X

k

)T

H

k

gk

(gk

)T

H

kk

1

k其中：X

k

X

k

1

X

kgk

f

(

X

k

1

)

f

(

X

k

)

BFGS

DFPSk

1

H

k

1f

(

X

k

1

)Matlab求解无约束规划问题类

型模

型基本函数名一元函数极小Min

F(x)s.t.x1<x<x2x=fminbnd(‘F’,x1,x2)多元无约束极小Min

F(X)X=fminunc(‘F’,X0)X=fminsearch(‘F’,X0)完整的调用格式：[X,FVAL,EXITFLAG,OUTPUT]

=fminbnd(FUN,x1,x2,OPTIONS,P1,P2,...)[X,FVAL,EXITFLAG,OUTPUT,GRAD,HESSIAN]

=fminunc(FUN,X0,OPTIONS,P1,P2,...)[X,FVAL,EXITFLAG,OUTPUT]

=fminsearch

(FUN,

X0,

OPTIONS,

P1,P2,...)输出变量含义变量描

述x由优化函数求得的值.若exitflag>0,则x为解;否则,x不是最终解,它只是迭代终止时优化过程的值fval解x处的目标函数值exitflag描述退出条件:exitflag>0,表目标函数收敛于解x处exitflag=0,表已达到函数计算或迭代的最大次数exitflag<0,表目标函数不收敛output包含优化结果信息的输出结构.Iterations:

迭代次数Algorithm:所采用的算法FuncCount:计算函数值次数AlgorithmsLarge-Scale

Optimization——

By

default

fminuncchooses

the

large-scale

algorithm

if

the

usersupplies

the

gradient

in

fun

(and

GradObj

is

'on'in

options).

This

algorithm

is

a

subspace

trustregion

method

and

is

based

on

the

interior-reflective

Newton

method.Medium-Scale

Optimization——fminunc

withoptions.LargeScale

set

to

'off'

uses

the

BFGSQuasi-Newton

method

with

a

mixed

quadratic

andcubic

line

search

procedure.

This

quasi-Newtonmethod

uses

the

BFGS

formula

for

updating

theapproximation

of

the

Hessian

matrix.

The

DFPformula,

which

approximates

the

inverse

Hessianmatrix is

selected

by

setting

optionsDisplay:显示水平.取值为’off’时,不显示输出;取值为’iter’时,显示每次迭代的信息;取值为’final’时,显示最终结果.默认值为’final’.MaxFunEvals:允许进行函数计算的最大次数,取值为正整数.MaxIter:

允许进行迭代的最大次数,取值为正整数.控制参数options可以通过函数optimset创建或修改.格式为options=optimset(‘param1’,value1,’param2’,value2,...)创建一个名称为options的优化选项参数,其中指定的参数具有指定值,所有未指定的参数取默认值.例：opts=optimset(‘Display’,’iter’,’TolFun’,1e-8)该语句创建一个称为opts的优化选项结构,其中显示参数设为’iter’,TolFun参数设为1e-8.options中常用的几个参数的名称、含义、取值如下:用Matlab解无约束优化问题1.一元函数无约束优化问题:常用格式如下：min

f（x）x1

x

x2x=

fminbnd(fun,x1,x2)x=

fminbnd

(fun,x1,x2

,

options)(3)[x，fval]=

fminbnd(...)(4)[x，fval，exitflag]=

fminbnd(...)(5)[x，fval，exitflag，output]=

fminbnd(...)函数fminbnd的算法基于0.618算法和二次插值法，它要求目标函数必须是连续函数，并可能只给出局部最优解。例

1

求 f

=

2

e

x

sin

x

在

0<x<8

中的最小值与最大值主程序为test.m:f='2*exp(-x).*sin(x)';f1='-2*exp(-x).*sin(x)';fplot(f,[0,8],'*-');

%作图语句hold

on,fplot(f1,[0,8],

'^-');

%作图语句

[xmin,ymin]=fminbnd

(f,0,8)[xmax,ymax]=fminbnd

(f1,

0,8)运行结果：xmin

=

3.9270xmax

=

0.7854ymin

=

-0.0279ymax

=

0.64482、多元函数无约束优化问题标准型为：min

F(X)命令格式为:x=fminunc(fun,X0

)；或x=fminsearch(fun,X0

)x=

fminunc(fun,X0

,

options)；或x=fminsearch(fun,X0

,options)(3)[x,fval]=fminunc(...)；或[x,fval]=fminsearch(...)(4)[x,fval,exitflag]=fminunc(...)；或[x,fval,exitflag]=

fminsearch(5)[x,fval,exitflag,output]=

fminunc(...)；或[x,fval,exitflag,output]=fminsearch(...)fminsearch是用单纯形法寻优.fminunc的算法见以下几点说明：fminunc为无约束优化提供了大型优化和中型优化算法。由

options中的参数LargeScale控制：LargeScale=’on’(默认值),使用大型算法(信赖域算法)LargeScale=’off’(默认值),使用中型算法fminunc为中型优化算法的搜索方向提供了4种算法，由

options中的参数HessUpdate控制a：HessUpdate=’bfgs’(默认值),拟牛顿法的BFGS公式；HessUpdate=’dfp’，拟牛顿法的DFP公式；

HessUpdate=’steepdesc’，最速下降法fminunc为中型优化算法的步长一维搜索提供了两种算法，由options中参数LineSearchType控制：

LineSearchType=’quadcubic’(缺省值)混合的二次和三次多项式插值；LineSearchType=’cubicpoly’，三次多项式插使用fminunc和fminsearch可能会得到局部最优解.例3

minf(x)=(4x1

+2x2

+4x1x2+2x2+1)*exp(x1)2

21、编写M-文件fun1.m:function

f=fun1

(x)f

=

exp(x(1))*(4*x(1)^2+2*x(2)^2+4*x(1)*x(2)+2*x(2)+1);2、输入M文件test.m如下:x0

=[-1,

1];x=fminunc(‘fun1’,x0);y=fun1(x)3、运行结果:-1.0000x=

0.5000y

=1.3029e-10例4 Rosenbrock函数

f(x1，x2)=100(x2-x12)2+(1-x1)2

的最优解(极小)为x*=(1,1)，极小值为f*=0.试用不同算法(搜索方向和步长搜索)求数值最优解.

初值选为x0=(-1.2

,

2).为获得直观认识,先画出Rosenbrock函数的三维图形,输入以下命令:[x,y]=meshgrid(-2:0.1:2,-1:0.1:3);z=100*(y-x.^2).^2+(1-x).^2;mesh(x,y,z)画出Rosenbrock

函数的等高线图,输入命令：contour(x,y,z,20)hold

onplot(-1.2,2,'

o

');text(-1.2,2,'start

point')plot(1,1,'o')text(1,1,'solution')3.用fminsearch函数求解输入命令:f='100*(x(2)-x(1)^2)^2+(1-x(1))^2';[x,fval,exitflag,output]=fminsearch(f,

[-1.2

2])运行结果:x

=1.0000

1.0000fval

=1.9151e-010exitflag

=

1output

=iterations:

108funcCount:

202algorithm:

'Nelder-Mead

simplex

direct

search'使用fminunc函数见Matlab演示！问题：增加生产、发展经济所依靠的主要因素有增加投资、增加劳动力以及技术革新等，试建立模型讨论究国民经济产值与这些因素的数量关系。简化：由于技术水平不像资金、劳动力那样容易定量化,初步认为技术水平不变,只讨论产值和资金、劳动力之间的关系(在技术发展不快的时期,该简化是有意义的.)分析及建模：用Q,K,L分别表示产值、资金和劳动力，该模型希望建立合适的函数Q(K,L)来描述三个量间的关系。假设1)产值依赖于每个劳动力的投资强度，并且与劳动力数示例:经济增长模型——最小二乘问题量成正比，即有(G表示某一函数)：Q(

K

,

L)

G(

y)L,

y

LK假设2)

产值随资金、劳动力的增加而增加，但增长得越来越慢.

据此选择G函数的形式为G(

y)

ayb

,

a

0,

0

b

1简化得:Q(

K

,

L)

aK

b

L1b——这是经济学中著名的柯布-道格拉斯(Cobb-Douglas)生产函数，其更一般的形式为:Q(

K

,

L)

aK

b

Lc

,

0

b,

c

1其中a,b,c由统计数据确定——最小二乘问题!现有统计数据(Qi,Ki,Li),i＝1,2,...,n.代入上式分别可得n个方程:b

ci

i

i

iQ

aK

L

0f

(a,

b,

c)i

1,2,,

n由于统计特性，一般来说，这些等式不可能完全成立，记:b

ci

i

i

iQ

aK

L最小二乘问题就是使得该误差项尽可能小，于是极小化：b

ci

i

iQ

aK

L

ni

1i2ni

12

mini

1,2,,

n这是非线性最小二乘问题！ni

1Q(

K

,

L)

aK

b

Lc

,

0

b,

c

1对该式两端取对数，可得ln

Q

ln

a

b

ln

K

c

ln

L

a

b

ln

K

c

ln

L这样对统计数据(Qi,Ki,Li),i＝1,2,...,n.代入上式分别可得n个线性方程:fi

(a,

b,

c)

lnQi

a

b

ln

Ki

c

ln

Li

0i

1,2,,

n使得以上等式误差项尽可能小，于是极小化：i

i

i2

c

ln

Lln

Q

a

b

ln

K

min这样就转化为线性最小二乘问题！——这更容易求解模型验证美国马萨诸赛州1890－1926年统计数据(以1899年归1)tQKLtQKLtQKL18900.720.950.7819031.301.221.2219162.093.611.8618910.780.960.8119041.301.271.1719171.964.101.9318920.840.990.8519051.421.371.3019182.204.361.9618930.730.960.7719061.501.441.3919192.124.771.9518940.720.930.7219071.521.531.4719202.164.751.9018950.830.860.8419081.461.571.3119212.084.541