gradientmethod.cpp

线性搜索方法。用梯度法求单峰极值。适合初学者参考

#include <iostream.h>
#include <math.h>
#include "LineSearch.h"
#define PRECI2 0.000001  // [ 精度要求 ]

double func( double *x );  // [ 目标函数的原型说明 ]


static double *start = new double[2];  // [ 初始点 ]
static int n;                          // [ 变量个数 ]  
static double *d = new double[n];      // [ 负最速下降方向 ]

// [ 用理查德外推法计算偏导数 ]
void ComputeGrad( double(*pf)(double *x),int num,double *point,double *grad )
{
	int i;  // [ 循环变量 ]
	double h = 1e-3;
	double *temp = new double[num];

	for( i=0; i<num; i++ )
	{
		temp[i] = point[i];
	}

	for( i=0; i<num; i++ )
	{
		temp[i] += 0.5*h;
		grad[i] = 4*pf(temp)/(3*h);
		temp[i] -= h;
		grad[i] -= 4*pf(temp)/(3*h);
		temp[i] += 3*h/2;
		grad[i] -= pf(temp)/(6*h);
	    temp[i] -= (2*h);
		grad[i] += pf(temp)/(6*h);
		temp[i] = point[i]; 
	}

	delete []temp;
}

// [ 求一向量的二范数 ]
double Norm( double *vector,int num )
{
	double norm = 0; 
	int i;   // [ 循环变量 ]

	for( i=0; i<num; i++ )
	{
		norm += vector[i]*vector[i];
	}

	norm  = sqrt(norm);
	return norm;
}

// [ 关于步长因子的函数 ]
double pfstep( double step )
{
	double *steppoint = new double[n];  // [ 跨步后的点 ]
	double value;   // [ 函数值 ]
	int i;   // [ 循环变量 ]

	for( i=0; i<n; i++ )
	{
		steppoint[i] = start[i] - step*d[i];
	}
	value = func(steppoint);
	delete []steppoint;
	return value;
}

//[ 梯度法求函数极小值 ]
double GradientMethod( double(*pf)(double *x),int num,double *start )
{
//	static double *d = new double[n];  // [ 最速下降方向 ]
	double step;                // [ 最佳步长 ]
	int i;   // [ 循环变量 ]

	ComputeGrad( pf,num,start,d); // [ 计算最速下降方向 ]

	while( Norm(d,num) > PRECI2 )
	{

		// [ 一维搜索求步长因子 ]
		step = LineSearch( pfstep,1 );

		// [ 计算新的迭代点 ]
		for( i=0; i<n; i++ )
		{
			start[i] = start[i] - step*d[i];
		}

		// [ 计算新的下降方向 ]
		ComputeGrad( pf,num,start,d);
	}

	// [ 返回极小值 ]
	return pf(start);
}

// [ 目标函数 ]
double func( double *x )
{
	return ( (x[0]-1)*(x[0]-1)+100*(x[1]-x[0]*x[0])*(x[1]-x[0]*x[0]));
//	return (100*(x[0]-1)*(x[0]-1) + (x[1]-1)*(x[1]-1));
//	return ( -2*x[0]+101*x[0]*x[0]-200*x[0]*x[1]+100*x[1]*x[1] );
}
// [ 主函数 ]
void main(void)
{
//	static double *start = new double[2];
	double min;
//	static int n;
	
	start[0] = -2;
	start[1] = 4;
//	start[0] = 0;
//	start[1] = 0;
    n = 2;

	min = GradientMethod( func,n,start );

	cout<<"极小点为：("<<start[0]<<","<<start[1]<<")"<<'\n';
	cout<<"极小值为："<<min<<'\n';

//	delete []start;
//	delete []d;
}