实验一工业生产的产量问题.docx

实验七　生产计划中的产量问题
【实验目的】
1．介绍无约束最优化方法的一些基本概念。
2．了解几种常见的无约束优化问题的求解方法，如迭代算法、最速下降法（梯度法）、牛顿法（Newton）、拟牛顿法。
3．学习掌握用MATLAB优化工具箱中的命令来求解无约束优化问题。
【实验内容】
某公司生产一种产品有甲、乙两个品牌，试讨论产销平衡下的最大利润。所谓产销平衡指公司的产量等于市场上的销量。利润既取决于销量和单件价格，也依赖于产量和单件成本。按照市场规律，甲种品牌的价格固然会随其销量的增长而降低；同时乙品牌销量的增长也会使甲的价格有稍微下降，根据实际情况，可以确定价格与销量成线性关系，即
　　　　　　　　＝300－2.35－0.09
乙的价格遵循同样的规律，有
＝480－0.14－2.98
甲品牌的成本会随着其产量的增长而降低，按实际情况可假设为负指数关系，即有
　　　　　　　　＝38＋116
乙品牌的成本遵循同样的规律，有
　　　　　　　　＝94＋145
试确定甲、乙两种品牌的产量，使公司获得的总利润最大。
【实验准备】
无约束最优化方法是指在没有约束条件限制下，求多变量实值函数极值的方法。无约束最优化问题的数学表达式为
　　　　　　　　，　＝（，，…，）　　　　　　（1）
一般为非线性函数，是维实变量，实际上这是一个多元函数无条件极值问题。由于一个求极大值问题，可以添加负号的方式转化为求极小值问题，因此通常只讨论求极小值问题。应该注意的是，极值问题的解，即极值点，都是局部最优解，全局最优解只能从局部最优解的比较中得到。
如何求解无约束最优化问题（1）的最优解呢？一般是采用迭代方法，即先选择一个初始点，再寻找该点处的下降方向，我们称为搜索方向。在该方向上求极小点，得到一个新的点。这个新的点要优于原来的点，即新点处的目标函数值小于原来点处的目标函数值。然后在新点处再寻找下降方向和该方向上求极小点，……，如此下去，最终得到最优点。
因此，求解无约束最优化问题需解决两个问题：一是在某些方向上的一维极小点，我们也称为一维搜索；另一个是寻找某些点处的下降方向，这是无约束最优化方法中最重要的一个问题。我们先了解第一个问题最常用的搜索方法。
1．求一维极小的二点二次插值方法
设是点处的一个搜索方向，要在该方向上寻优问题，转化为求一维函数
　　　　　　　　＝（＋）　　　　　　　　　　　　　　　 （2）
的求极值问题。
最常用的一维搜索方法是插值法，即用某些点的函数值（或导数值）构造插值函数，用插值函数的极小点来近似函数的极小点。这里介绍一种有效的插值方法，称为二点二次插值方法，即用二点处的函数值和一个点处的导数值构造二次函数，反复用二次函数的极小点来逼近函数的极小点。
算法1　二点二次插值法
（1）取初始点（如＝0），初始步长（如＝1）和步长缩减因子（0，1）（如＝0.1），置精度要求，并计算＝，＝
（2）若＜0，则置＝；否则置＝－
（3）计算＝＋和＝
（4）如果≤＋（－），则置＝2，转（3）
（5）计算＝－，＝，＝（）
（6）若≤，则停止计算（作为极小点）；否则置＝，＝，＝，＝，转（2）
2．最速下降法
前面介绍了一维搜索的二点二次插值方法，下面讨论如何选择搜索方向的问题，我们先来看看两个概念。
定义1　称维向量（，，…，）T为函数在处的梯度，记为▽
　　　　　　　　▽＝（，，…，）T　　　　　　　　　　　（3）
定义2　设是任意的单位向量，若极限存在，则称该极限为函数在处沿方向的一阶方向导数，简称为方向导数，记为，即
　　　　　　　　＝　　　　　　　　　　　　（4）
最速下降法的基本思想：选取一点作为初始点，计算该点的梯度▽，求该点处的最速下降方向，即令＝－▽，再沿方向前进，寻找该方向上的极小点，得到点，再计算▽，令＝－▽，沿方向前进，得到点，如此下去……具体算法如下：
算法2　最速下降法
（1）取初始点，置精度要求为，置＝1
（2）若||▽ ||＜，则停止计算（为无约束问题的最优解）；否则置
＝－▽
（3）一维搜索，求一维问题
　　　　　　　　＝（）
得，置＝。
（4）置＝＋1，转步骤（2）
在算法中，为精度要求，即当梯度接近于0时，我们就认为达到极小点，终止计算。这样做的目的是避免算法产生死循环，算法中的一维搜索可用算法1来求解。
最速下降法是一种最基本的算法，它在最优化方法中占有重要地位。其优点是工作量小，存储变量较少，初始点要求不高；缺点是收敛慢，效率不高，有时达不到最优解。下面介绍一种简单而直观的方法——Newton法。
3．牛顿法
将的二阶导数记作▽2＝（）是一个矩阵，称黑塞（Hessian）矩阵（简记）。当二阶偏导数连续时，是对称矩阵。
Newton法的基本思想是用一个二次函数去近似目标函数，然后精确地求出这个二次函数的极小点，并把这个极小点作为所示函数的极小点的近似值。
设是极小点的第次近似，将在点作二阶Taylor展开，略去高于二阶的项，用＝得到
＝≈＋▽＋▽2　 　　　 （5）
为了求使最小，只需将（5）式右端对求导，并令其为零，可得
　　　　　　　　▽＋▽2＝0　　　　　 　　　　　　　　 （6）
称为牛顿方程。当黑塞矩阵▽2正定时，可解得
　　　　　　　　＝－(▽2)-1▽　　　　　　　　　　　　　 （7）
称为牛顿方向。其迭代公式为
　　　　　　　　　　＝－(▽2)-1▽　　　　　　　　　　　　（8）
4．拟牛顿法
牛顿法的优点是在接近最优解时具有2阶收敛性，缺点是计算复杂，而且会出现病态、不正定等情况。为克服这些问题，同时又保持较快的收敛，利用第和第步得到的、、▽、▽根据拟牛顿条件构造一正定矩阵近似代替▽2 或近似代替(▽2)-1，得到的迭代算法称为拟牛顿法。
其中有两个常用的著名的迭代公式为BFGS公式和DFP公式，可参见其他最优化方法书籍或MATLAB帮助。
5．MATLAB求解无约束优化问题的主要命令
（1）MATLAB5.2及以下版本使用的命令

x = fmin( ' fun ' , x1 , x2 )　 求一元函数y = f( x )在[ x1 , x2 ]内的极小值
x = fmin( ' fun ' , x1 , x2 , options )　同上，参数options的定义由表1给出
[ x , options ] = fmin( ... )　同上，同时返回参数options的值
fmin函数采用黄金分割法和抛物线插值法，fun可直接用函数表达式表示，也可以是用M文件定义的函数名。
建立函数形式的M–文件格式如下，文件名fun.m（一般M函数文件名取函数名）
　　　　function f = fun( x )
f = f( x )

x = fmins( ' fun ' , x0 )　求多元函数（1）以x0为迭代初值的局部极小值
x = fmins( ' fun ' , x0 , options )　同上，参数options的定义由表1给出
[ x , options ] = fmins( ... )　同上，同时返回参数options的值
fmins函数采用Nelder–Meade单纯形搜索法，fun可直接用函数表达式表示，也可以是用M文件定义的函数名。
x = fminu( ' fun ' , x0 )　求多元函数（1）以x0为迭代初值的局部极小值
x = fminu( ' fun ' , x0 , options )　同上，参数options的定义由表1给出
[ x , options ] = fminu( ' fun ' , x0 , ... )　同上，同时返回参数options的值
fminu函数为无约束优化提供了三种算法，由options（6）控制，为步长一维搜索提供了二种算法，由options（7）控制。这里，fminu必须先用M–文件定义函数fun。
（2）MATLAB5.3以上版本使用的命令
fminbnd替代fmin求解一元函数极值，使用格式、搜索算法与之相同，[ x , Fval ] = fminbnd( ' fun ' , x0 ,... )同时返回解x处的函数值，而不是参数options的返回值
fminsearch替代fmins求解多元函数（1）极值，使用格式、搜索算法与之相同
fminunc替代fminu求解多元函数（1）极值，使用格式、搜索算法与之相同
有关上述命令的详细信息和使用方式可在帮助文件中了解，或在命令框里输入help “命令名”查阅。在MATLAB5.3以上的各种最新版本中fmin、fmins和fminu命令仍然有用，但在MATLAB的未来版本将可能删除这些命令。
（3）命令中参数options的有关定义
在大多数MATLAB优化命令函数中有一个控制参数options，它是一个有18个分量的向量，包含了在优化程序中要用到的参数，以便在计算最优值时控制精度要求、输出形式、搜索算法、迭代次数、步长等等是。在对优化命令函数的第一次调用时，向量options会自动使用缺省值；当然在调用前，可以对options的某些分量进行赋值，以达到控制要求。
表1　参数options各分量说明
序号
功能定义
分量说明
缺省值含义
1
输出形式
options(1)=1，有中间结果输出
options(1)=－1，给出警告信息
0，无中间结果输出
2
解的精度
options(2)设置的精度
≤－4
3
函数的精度
options(3)设置的精度
≤－4
4
函数的精度
options(4)设置的精度，由命令constr、minmax使用终止判据
≤－7
5
主要算法
options(5)=1，高斯–牛顿法
0，LM方法
6
搜索方向算法
options(6)=1，拟牛顿法DFP公式
options(6)=2，最速下降法
0，拟牛顿法的BFGS公式
7
步长一维搜索算法
options(7)=1，三次多项式插值
0，混合的二次和三次多项式
8
函数值输出
options(8)输出极值点处的函数值

9
梯度检查
options(9)=1，最初几步，梯度将与有限微分计算的结果比较
0，不作检查
10
函数计算次数
option(10)输出函数计算次数

11
梯度计算次数
option(11)输出梯度计算次数

12
约束梯度计算次数
option(12)输出约束梯度计算次数

13
等式约束的个数
option(13)确定等式约束数目
0，等式约束个数为0
14
最大迭代次数
option(14)输入迭代的最大次数
100，为自变量个数
200（fmins）
500（fmin）
15
目标优化
由attgoal命令使用，option(15)输入须达目的的目标个数
0
16
差分步长最小值
用差分计算梯度时步长的下限
－8
17
差分步长最大值
用差分计算梯度时步长的上限
0.1
18
步长
步长参数，一般取1或更小

【实验方法与步骤】
1．MATLAB优化工具箱中解无约束问题的命令和参数options的基本用法
下面用例题来予以说明
例1　求解 ＝－5，1≤≤2
输入命令：
>> x=fmin( 'exp(x) - 5*x' , 1 , 2 ) , f=exp( x ) - 5*x
x =
1.6094
f =
-3.0472
同时，该问题可以用fminbnd求解，得到的解答一样
输入命令：
>> [x , fval]=fminbnd( 'exp( x ) - 5*x' , 1 , 2 )
x =
1.6094
fval =
-3.0472
例2　用拟牛顿法的DFP公式求解＝＋－，初值为（1，2），要求精度为－7，给出函数计算次数以及函数值
首先建立M–文件，文件名取函数名fun.m
function f = fun( x )
f = 4 * x( 1 )^2 + x( 2 )^2 - x( 1 )^3 * x( 2 )
输入命令：
>> x0=[1, 2];
>> opt( 2 )=1e-6;　　%设置自变量要求的精度
>> opt( 3 )=1e-6;　　%设置函数所要求的精度
>> opt( 6 )=1;　　%搜索算法用拟牛顿法的DFP公式求解
>> [x, opt]=fminu( ' fun ', x0, opt);　　%返回参数options的值给向量opt
可得到结果
x =
1.0e-008 *
-0.3493 0.7566
参数options的值返回给向量opt，全部元素略。
>> n=opt( 10 );　　%函数计算次数赋值给n
n =
35
>> y=opt( 8 );　　%在极值点处函数值赋值给y
y =
1.0606e-016
2．引例问题的分析与求解
由题设可知，甲品牌产品单件获利为－，乙品牌产品单件获利为－，由产销平衡原理，所有产品的销量即为产量，则甲、乙两种产品总获利为
　　　　　　　　＝（－）＋（－）
容易看出，原问题实际上转化为求二元函数的极大值，为用MATLAB优化工具箱中的fminunc求解，需将其转化为求函数＝－的最小值。
为确定初始值，先忽略成本，并令价格中项的较小系数0.09和中项较小的系数0.14等于零（因为它们对价格的作用比较微弱，暂时可忽略不计），则确定初值问题转化为求
＝（300－2.35）＋（480－2.98）
的极值，很容易可以求得＝＝63.83，＝＝80.54，我们用它作为原问题的初始解。
首先建立M–文件，文件名取函数名fun1.m
function y = fun1( x )
p1 = 300 - 2.35 * x( 1 ) - 0.09 * x( 2 );
q1 = 38 * exp( - 0.023 * x( 1 ) ) + 116;
p2 = 480 - 0.14 * x( 1 ) - 2.98 * x( 2 );
q2 = 94 * exp( - 0.018 * x( 2 ) )+145;
y = - ( p1 - q1 ) * x( 1 ) - ( p2 - q2 ) * x( 2 )
　　输入命令：
>> x0=[63.83; 80.54];
>> [x, fval]=fminunc( ' fun1 ', x0);
可得到结果：
x =
35.8482
54.7380
fval =
-1.0015e+004
甲种品牌产量为35.8482，乙种品牌产量为54.7380，最大利润和为1.0015e+004。用fminu可得到同样结果，同时根据参数options可观察到函数采用拟牛顿的DFP公式并迭代计算了20次。
【练习与思考】
1．由于市场竞争的影响，电视机售价越高，销售量就会越低，它们之间关系为
　　　　　　　　　　＝，　，＞0
其中为最大需求量，为价格系数。另一方面，销售量越大，每台电视机成本就越低，
＝－，　，＞0，＞1
其中是只生产一台电视机的成本，为规模系数。应如何确定电视机售价才能获得最大的利润。