最优化课程讲义

最优化理论，第八章 无约束最优化的直接法
解析法 vs  直接法

产生背景
解析法，需要计算目标函数的一阶导数，甚至二阶导数。
目标函数的解析表达式比较复杂，或者难以用明显的解析式表达出来时，其导数难求或者无法求出。
？
给出一些只涉及目标函数值的计算、不涉及目标函数求导的求解方法。
仅仅利用目标函数值的信息直接建立搜索求解的方法通常称为直接搜索法（direct search method），简称直接法。
尽管算法的效果有时候不理想，但是构思直观、使用方便、效果稳定，实际工作者很愿意采用。

8.1坐标轮换法
8.2模式搜索法
8.3旋转方向法
8.4 Powell法
8.5单纯形调优法

讲述要点：基本思想、操作步骤、编程实现、方法评价

总体比较：

中文名，英文名 提出作者（时间）	一句话概括
坐标轮换法 Univariate search technique D’esopo （1995）	每一次迭代，以一个变量的坐标轴方向作为搜索方向，将N维的优化问题转化为一维搜索问题。
模式搜索法 Pattern search method Hooks 和 Jeeves（1961）	每一次迭代，交替进行轴向移动和模式移动。 轴向移动：探测下降的有利方向 模式移动：沿着有利方向加速移动
旋转方向法 Rotating direction method Rosenbrock（1960）	每一次迭代，采用变步长的轴向移动，然后利用轴向的旋转产生一组新的方向，作为下一次迭代的轴向。

8.1坐标轮换法

做一维搜索，究竟应该把哪个坐标轴方向放在前面，哪个坐标轴方向放在后面？
这会对算法的效率产生较大的影响！
实际问题中，可按照各个因素（变量）对实验结果影响的大小，依次从前往后排列坐标轴方向。

（分析各个因子对实验结果影响的大小，方法有层次分析法、专家打分法等。 大家懂的……）

算法描述

方法评价：
简单，容易实现。
当维数增加时，效率明显下降。
收敛慢，以振荡方式逼近最优点。
受目标函数的形态影响很大

8.2模式搜索法

首先，探测有利的下降方向：轴向移动（axis direction move）

然后，沿着有利方向加速移动：模式移动（pattern move）

算法描述

                                                      输入                               输出
说明                             初始值  收缩因子 要求精度     迭代次数      自变量X             函数值
书上的                          0，0     0.25       0.1             2               0.375 , 0.125     -1.016
改变初始值                    1，1     0.25       0.1             2               1.375 , 0.750     -2.703（好点）
改变精度                       0，0     0.25       0. 001        8               0.498 , 0.248     -1.308（不好）

改变初始值和收缩因子     1，1     0.5        0.1             3                1.875 , 1.000     -2.984（好点）

三个都改                       1，1     0.5        0.001          9               1.998 , 1.000     -3.000（很好）

算法评价：
可以体会到，选取初值是个大问题，确定收缩因子也比较棘手。
通过探测移动和模式移动，逐步逼近最优点，需要计算很多次函数值。（计算机不怕麻烦的~）

8.3旋转方向法

变步长轴向移动
轴向旋转，产生新的轴向
“旋转方向法，Rotating direction method ，是Rosenbrock在1960年提出来的。”
图书馆：没找到“旋转方向法”
陈开周，最优化计算方法，西北电讯工程学院出版社，1985.10   P253  Rosenbrock坐标轮换法，于1960年提出。
陈开周：每次迭代分为两部分，
(1)沿n个正交方向进行试探

(2)决定n个新的正交方向以改善原搜索方向，每次迭代的试探方向组是不同的。每一次迭代，采用变步长的轴向移动，然后利用轴向的旋转产生一组新的方向，作为下一次迭代的轴向。

算法描述

代码1：一维搜索之0.618法

P140
用0.618法求解    P136
#include<stdio.h>
#include<math.h>

//解析式原型
float f(float x)
{
	return (x*x*x-2*x+1);
}

//0.618法
float f_6_2()
{
	float a=-2;//低区间
	float b=2;//高区间
	float epsilon=0.001 ;//精度
	float tao=0.618;//tao的值

	float lambda=a+(1-tao)*(b-a);
	float miu=a+tao*(b-a);

	while(fabs(lambda-miu)>=epsilon)
	{
		float f1=f(lambda);//低者的函数值
		float f2=f(miu);//高者的函数值
		if(f1<f2)//向左搜索
		{
			b=miu;
			miu=lambda;
			lambda=a+(1-tao)*(b-a);
		}
		else//向右搜索
		{
			a=lambda;
			lambda=miu;
			miu=a+tao*(b-a);
		}
	}
	float t=(lambda+miu)/2.0;
	return t;
}

int main()
{
	printf("%f\n",f_6_2());
	return 1; 
}

运行结果：0.817043
参考P137，该结果比较好

代码2：模式搜索法

// 模式搜索法.cpp : 定义控制台应用程序的入口点。
//

//尝试探测每一个方向；若下降了，采取下降的建议；否则不采取建议。

#include "stdafx.h"
#include <iostream>
using namespace std;

double f(double x1,double x2)		//问题的解析式
{
//	return ((x1-1)*(x1-1)+5*(x1*x1-x2)*(x1*x1-x2));
	return (x1*x1+x2*x2-3*x1-x1*x2);
}

int _tmain(int argc, _TCHAR* argv[])
{
	double x1=1;			//初始值
	double x2=1;			//初始值
	double t1=x1;			//探测的初始值，变量x1的探头
	double t2=x2;			//探测的初始值，变量x2的探头
	double alpha=0.5;		//收缩因子？？
	double epsilon=0.001;	//精度
	double temp_x1;
	double temp_x2;

	int times=1;			//迭代次数
	do 
	{
		printf("第%2d 轮计算：",times);
		//第一个方向
		if (f(t1+alpha,t2)<f(t1,t2))
		{
			t1=t1+alpha;
		}
		if (f(t1-alpha,t2)<f(t1,t2))
		{
			t1=t1-alpha;
		}
		//第二个方向
		if (f(t1,t2+alpha)<f(t1,t2))
		{
			t2=t2+alpha;
		}
		if (f(t1,t2-alpha)<f(t1,t2))
		{
			t2=t2-alpha;
		}
		//模式移动
		if (f(t1,t2)<f(x1,x2))
		{
			temp_x1=x1;
			temp_x2=x2;
			x1=t1;
			x2=t2;
			t1=x1+x1-temp_x1;
			t2=x2+x2-temp_x2;
		}
		//检测上述的探测是否成功。若成功，继续；否则，回滚
		if ((x1!=t1||x2!=t2))
		{
			t1=x1;
			t2=x2;
		}
		else
		{
			break;
		}
		printf("t=( %5.3f , %5.3f )\n",t1,t2);
		//缩短步长
		alpha=alpha/2;			
		times++;
	} while (alpha>=epsilon);

	printf("结果值：%5.3f\n",f(t1,t2));

	system("pause");
	return 0;
}

代码3：旋转方向法

// 旋转方向法.cpp : 定义控制台应用程序的入口点。
//

#include "stdafx.h"
#include <iostream>
#include <math.h>
using namespace std;


double f(double x1,double x2)
{
	return (x1*x1+x2*x2-3*x1-x1*x2);
}
int main(void)
{
	double x1=0;
	double x2=0;
	double alpha=0.5;	//收缩因子
	double beta=3;		//放大因子
	double delta1=1;	//x1的步长，初始为1
	double delta2=1;	//x2的步长，初始为1
	double epsilon=0.0001;

	double y1;	//x1的参考点
	double y2;	//x2的参考点
	double z1;	//暂时保存z1的值
	double z2;	//暂时保存z2的值


	y1=x1;
	y2=x2;
	z1=x1;
	z2=x2;
	for(int times=1;times<=10;times++) 
	{

		cout<<endl<<"第"<<times<<"次迭代:"<<endl;
		//第一个方向
		if (f(y1+delta1,y2)<f(y1,y2))
		{
			y1=y1+delta1;
			delta1=delta1*beta;//放大
			cout<<"第一个方向，放大"<<beta<<"倍，delta1="<<delta1<<endl;
		}
		else
		{
			delta1=delta1*(-alpha);//缩小
			cout<<"第一个方向，缩小"<<(-alpha)<<"倍，delta1="<<delta1<<endl;
		}
		//第二个方向
		if (f(y1,y2+delta2)<f(y1,y2))
		{
			y2=y2+delta2;
			delta2=delta2*beta;//放大
			cout<<"第二个方向，放大"<<beta<<"倍，delta2="<<delta2<<endl;
		}
		else
		{
			delta2=delta2*(-alpha);//缩小
			cout<<"第二个方向，缩小"<<(-alpha)<<"倍，delta2="<<delta2<<endl;
		}
		if (f(y1,y2)<f(z1,z2))
		{
			z1=y1;
			z2=y2;
			cout<<"下降方向，更新解的值"<<endl;
		}
	} 
	cout<<endl<<"最终结果 x= "<<z1<<" , "<<z2<<endl;
	system("pause"); 
	return 0;
}
结果：
 
第1次迭代:
第一个方向，放大3倍，delta1=3
第二个方向，缩小-0.5倍，delta2=-0.5
下降方向，更新解的值
 
第2次迭代:
第一个方向，缩小-0.5倍，delta1=-1.5
第二个方向，缩小-0.5倍，delta2=0.25
 
第3次迭代:
第一个方向，缩小-0.5倍，delta1=0.75
第二个方向，放大3倍，delta2=0.75
下降方向，更新解的值
 
第4次迭代:
第一个方向，放大3倍，delta1=2.25
第二个方向，放大3倍，delta2=2.25
下降方向，更新解的值
 
第5次迭代:
第一个方向，缩小-0.5倍，delta1=-1.125
第二个方向，缩小-0.5倍，delta2=-1.125
 
第6次迭代:
第一个方向，缩小-0.5倍，delta1=0.5625
第二个方向，缩小-0.5倍，delta2=0.5625
 
第7次迭代:
第一个方向，缩小-0.5倍，delta1=-0.28125
第二个方向，缩小-0.5倍，delta2=-0.28125
 
第8次迭代:
第一个方向，缩小-0.5倍，delta1=0.140625
第二个方向，缩小-0.5倍，delta2=0.140625
 
第9次迭代:
第一个方向，放大3倍，delta1=0.421875
第二个方向，缩小-0.5倍，delta2=-0.0703125
下降方向，更新解的值
 
第10次迭代:
第一个方向，缩小-0.5倍，delta1=-0.210938
第二个方向，放大3倍，delta2=-0.210938
下降方向，更新解的值
 
最终结果 x= 1.89063 ,0.929688
 
 
循环次数越多，值越准确。

收获：

好书，启发思维；差书，混淆思路。