ch21-半无限问题

第21章 半无限问题

21.1 基本数学原理

半无限约束问题的数学模型如下所示：

),(...

),(0

),(0

)(0

)()(min 2211≤≤≤≤≤=≤=≤n n x w x K w x K w x K ub

x lb beq x Aeq b

Ax x ceq x c x f 其中，x,b,beq,lb 和ub 为向量，A 和Aeq 为矩阵，c (x ),ceq (x ) 和),(i i w x K 为返回向量的函数。f (x ) 为返回标量的函数。f(x ),c (x )和ceq (x ) 可以是非线性函数。向量（或矩阵）?)1(0),(11i w x K 是否这里的≤即为半无限约束，它是x 和其他变量n w w w ,...,,21的连续函数。n w w w ,...,,21变量的最大长度为2。

MATLAB 优化工具箱采用二次、三次混合插值法结合序列二次规划法(SQP)进行问题求解。

21.2 有关函数介绍

用fseminf 函数求半无限约束多变量非线性函数的最小值。其调用格式为： ●

fseminf 函数求几个变量半无限约束标量函数的最小值，初值给定。目标是使f(x)最小化。因为不可能计算所有的可能值，所以必须选择一个区域，在它上面计算样本数据集。 ●

x=fseminf(fun,x0,ntheta,seminfcon) 初值为x0,求约束条件为ntheta,半无限约束为seminfcon 的fun 函数的最小值。 ●

x=fseminf(fun,x0,ntheta,seminfcon,A,b)该函数试图满足线性不等式A*x<=b ●

x=fseminf(fun,x0,ntheta,seminfcon,A,b,Aeq,beq)在上面的基础上添加线性等式Aeq*x=beq. 当没有不等式存在时，设置A=[ ]、b=[ ] 。 ●

x=fseminf(fun,x0,ntheta,seminfcon,A,b,Aeq,beq,lb,ub)定义设计变量x 的一系列下界lb,和上界ub,使得总有Ib

x=fseminf(fun,x0,ntheta,A,b,Aeq,beq,Ib,ub,options)用options 结构指定的优化参数进行最小化。 ● x=fseminf(fun,x0,ntheta,seminfcon,A,b,Aeq,beq,Ib,ub,options,P1，P2…)将问题参数P1,P2

等直接传递给函数fun 和seminfcon.若变量A,b,Aeq,beq,Ib,Ub,和options 不需要，则将它们设置为空矩阵。

●[x,fvaI]=fseminf(…)返回解x处的目标函数值。

●[x,fvaI,exitfiag]=fseminf(…)返回描述退出条件的exitfiag参数。

●[x,fvaI,同esitfIag,output]=fseminf(…)返回包含优化信息的输出参数output。

●[x,fvaI,exitfIag,output,Iambda]=fseminf(…)返回包含解x处拉格朗日乘子的Iambda参数。

各调用格式中，ntheta参数为半无限约束的个数。options参数为优化参数选项，有以下一些属性：

DerivativeCheck——比较用户提供的导数（梯度）和有限差分导数。

Diagnostics——打印待最小化或待求解函数的诊断信息。

DiffMaxChange——变量有限差分的最大变化。

DiffMinChange——变量有限差分的最小变化。

Display——显示水平。设置为’off’时不显示输出; 设置为’iter’时显示每一步迭代过程的输出；设置为’final’时显示最终的迭代结果。

GradObj——用户定义的目标函数的梯度信息。对于大型算法，必须提供梯度；对于中型算法，梯度则是可选项。

MaxFunEvaIs——函数评价的最大次数。

MaxIter——函数迭代的允许最大次数。

ToICon——约束矛盾的终止容限。

ToTFun——函数值的终止容限，

ToIX——解x处的终止容限。

Seminfcon参数计算非线性不等式向量c，非线性等式向量ceq和ntheta半无限约束（向量或矩阵）K1,K2,…,Kntheta. seminfcon参数是一个包含函数名的字符串，该函数可以是M 文件、内部文件或MEX文件。例如，若seminfcon=’myinfcon’, 则M文件myinfcon.m有下面的形式：

function [c,ceq,K1,K2,…,Kntheta,s]=myinfcon(x,s) % 初始化样本区间

if isnan(S(1,1)),

S=…% S有ntheta行、2列

end

w1=…% 计算样本集

w1=…% 计算样本集

w2=…% 计算样本集

…

wntheta=…% 计算样本集

K1=…% x和w处的一阶半无限约束

K2=…% x和w处的二阶半无限约束

…

Kntheta=…% x和w处的最后一阶半无限约束

c=…% 计算x和w处的非线性不等式

ceq=…% 计算x和w处的非线性等式

S为建议的样本区间，它可能会被利用，也可能不会被利用。如果没有此约束存在，将为c和ceq返回空矩阵。

向量或矩阵K1,K2,…,Kntheta包含分别为样本集中的独立变量w1,w2,…,wntheta进行评价的半无限约束。两个列矩阵S包含w1,w2,…,wntheta值的建议的样本区间，它们用于评价K1,K2,…,Kntheta。S矩阵的第i行包含评价Ki的样本区间。当Ki为一向量时，只用S(i,j)（第2列可以都为零）。当Ki为一矩阵时,S(i,2)用于对Ki中的行进行取样,S(i,1)被用作Ki

中列的取样区间。第一次迭代，S 为空值，所以有些初始取样区间必须由seminfcon 参数来决定。

其他参数的意义可以参见表15-7和表15-8。

fseminf 函数使用二次和三次混合插值法估计半无限约束范围内的峰值。峰值用于组成一系列约束提供给SQP 法，就像在fmincon 函数中所做的一样。当约束的数目改变时，拉格朗日乘子被重新分配到新的约束集合中去。

建议计算取样区间时使用内插峰值和峰值之间的差值去估计插值点的数目。内插效果通过曲线外推和与曲线中其他点比较而被考虑。当峰值接近约束边界时，建议的取样区间便窄了。

使用fseminf 函数时需要注意下面的一些问题：

Seminfcon 参数中设置建议取样区间S 可能会被计算过程中的优化方法改变，因为

其他不是建议值的值可能会更有效或更稳健。所以，计算Ki(x1,wi)时，只要不引起局部极小值的明显变化，其中的wi 是允许变化的。

目标函数和约束函数，以及半无限约束函数必须为 x 和w 的连续函数。fseminf

函数可能只给出局部最优解。

当问题不可行时，fseminf 函数试图使最大的约束值最小化。

21.3 应用实例

【例21-1】一维问题。

假设有下面的函数：

232221)5.0()5.0()5.0()(-+-+-=x x x x f

式中

100

1 100

1 :1)sin()50(1000

1)cos()sin(),( 1)sin()50(1000

1)cos()sin(),( 21213322212222233121211111≤≤≤≤≤----=≤----

=w w w w x x w w x w x w w x K x x w w x w x w w x K 有下面的约束和 注意：这里半无限约束是一维向量——Ki 只依赖于wi,。因为约束条件必须写为0),(≤i i w x K 的形式，故需要对约束做如下计算。

首先，写一个M 文件,计算目标函数值。

function f=myfun2101(x,s)

f=sum((x-0.5).^2);

然后，编写M 文件,计算非线性等式约束和非线性不等式约束，以及半无限约束： function [c,ceq,K1,K2,s]=mycon2101(x,s)

% 初始化样本空间

if isnan(s(1,1)),

s=[0.2 0;0.2 0];

end

% 样本集

w1=1:s(1,1):100;

w2=1:s(2,1):100;

% 半无限约束

K1=sin(w1*x(1)).*cos(w1*x(2))-1/1000*(w1-50).^2-sin(w1*x(3))-x(3)-1; K2=sin(w2*x(2)).*cos(w2*x(1))-1/1000*(w2-50).^2-sin(w2*x(3))-x(3)-1; % 无约束

c=[];

ceq=[];

% 绘制半无限约束图

plot(w1,K1,'-',w2,K2,':'),title('Semi-infinite constraints')

drawnow

然后调用优化过程：

>> x0=[0.5;0.2;0.3];

>> [x,fval]=fseminf(@myfun2101,x0,2,@mycon2101)

经过8(?)次跌代以后，得到问题道的解: 解x及x处的函数值Optimization terminated successfully:

Search direction less than 2*options.TolX and

maximum constraint violation is less than options.TolCon

Active Constraints:

7

10

x =

0.6673

0.3013

0.4023

fval =

0.0770

半无限约束的最大值为

>> [c,ceq,K1,K2]=mycon2101(x,NaN);

>> max(K1)

ans =

-0.0017

>> max(K2)

ans =

-0.0845

生成半无限约束图，如图21-1所示。

该图演示了约束边界上两个函数如何达到峰值。

‘mycon2101.m 文件内部的绘图语句将使计算速度降低，剔除它可以加快速度。

【例21-2】 二维问题。

该问题的数学模型为

232221)2.0()2.0()2.0()(-+-+-=x x x x f % 书上0.2为0.5，但做不出来 式中

100

1 100

1 :5.1)sin()50(1000

1)cos()sin( ...)sin()50(1000

1)cos()sin(),( 21213322211223312122111≤≤≤≤≤----+----

=w w w w x x w w x w x w x x w w x w x w w x K 有下面的约束和 初始点为x=[0.2 0.2 0.2 ].

注意半无限约束是二维的，即为矩阵。

首先，写一个M 文件,计算目标函数值：

function f=myfun2102(x,s)

f=sum((x-0.2).^2);

第2步，为约束条件编写M 文件:

function [c,ceq,K1,s]=mycon2102(x,s)

% 初始化取样区间

if isnan(s(1,1)),

s=[2 2];

end

% 样本集

w1x=1:s(1,1):100;

w2y=1:s(1,2):100;

[wx,wy]=meshgrid(w1x,w2y);

%

% 半无限约束

K1=sin(wx*x(1)).*cos(wy*x(2))-1/1000*(wx-50).^2-sin(wx*x(3))-x(3)+... sin(wy*x(2)).*cos(wx*x(1))-1/1000*(wy-50).^2-sin(wy*x(3))-x(3)-1.5; %

% 无有限非线性约束

c=[];

ceq=[];

%

% 网络图

mesh(wx,wy,K1)

title('Semi-infinite constraints')

drawnow

%

% 也可用下面的4行作网络图，效果更好看

% m = surf(wx,wy,K1,'edgecolor','none','facecolor','interp');

% camlight headlight

% title('Semi-infinite constraint')

% drawnow

下一步调用优化过程M文件

>> x0=[0.25,0.25,0.25];

>> [x,fval]=fseminf(@myfun2102,x0,1,@mycon2102)

经过7次跌代以后，得到问题的解：

Optimization terminated successfully:

Magnitude of directional derivative in search direction

less than 2*options.TolFun and maximum constraint violation is less than options.TolCon

Active Constraints:

16

x =

0.2081 0.2066 0.1965

fval =

1.2126e-004

解x处的半无限约束的最大值为

>> [c,ceq,K1]=mycon2102(x,NaN);

>> max(max(K1))

ans =

-0.0260

生成网格图21-2所示。

图21-2 二维问题的半无限约束图