最速下降法求解无约束最优化问题
1 理论基础
已知问题模型为
()min n x R
f x ∈ 算法:
1) 选取初始点0x 与()00f f x =,初始化0k =
2) 验证迭代中止条件:k x eps ∇<且()k f x eps ∇<或lim k > 若是k x eps ∇<且()k f x eps ∇<,则输出解k x 及迭代次数k 若是lim k >,则输出error 信息
注:其中1010,lim 1000eps -==
3) 计算k x 点搜索方向:()k k d f x =-∇
计算k x 点迭代步长:(){}0
arg min k k k f x d ααα≥=+ 更新点列:1k k k k x x d α+=+,1k k =+
转第2步
2 MATLAB 程序
2.1 函数说明
文件名称:Opt_Steepest.m
function [xo, fo] = Opt_Steepest(f, grad, x0, TolX, TolFun, dist0, MaxIter)
% 用最速下降法求最优化解
% 输入:
% f——函数名
% grad——梯度函数
% x0——解的初值
% TolX——变量的误差阈值
% TolFun——函数的误差阈值
% dist0——初始步长
% MaxIter——最大迭代次数
% 输出:
% xo——最小值的点
% fo——最小的函数值
%% 判断输入的变量数,设定一些变量为默认值if nargin < 7
MaxIter = 1000; %最大迭代次数默认为100 end
if nargin < 6
dist0 = 10; %初始步长默认为10
end
if nargin < 5
TolFun = 1e-10; %函数值误差为1e-8
end
if nargin < 4
TolX = 1e-5; %自变量距离误差
end
%% 第一步,求解的初值的函数值
x = x0;
fx0 = feval(f, x0);
fx = fx0;
dist = dist0;
kmax1 = 25; % 线性搜索法确定步长的最大搜索次数
warning = 0;
%% 迭代计算求最优解
for k = 1 : MaxIter
g = feval(grad, x);
g = g/norm(g); % 求在x处的梯度方向
%线性搜索方法确定步长
dist = dist*2; % 令步长为原步长的二倍
fx1 = feval(f, x-dist*2*g);
for k1 = 1 : kmax1
fx2 = fx1;
fx1 = feval(f, x-dist*g);
if fx0 > fx1+TolFun && fx1 < fx2 - TolFun % fx0 > fx1 < fx2,
den = 4*fx1 - 2*fx0 - 2*fx2;
num = den - fx0 + fx2; % 二次逼近法
dist = dist*num/den;
x = x - dist*g;
fx = feval(f,x); % 确定下一点
break;
else
dist = dist/2;
end
end
if k1 >= kmax1
warning = warning + 1; %无法确定最优步长
else
warning = 0;
end
if warning >= 2 || (norm(x - x0) < TolX && abs(fx - fx0) < TolFun) break;
end
x0 = x;
fx0 = fx;
end
xo = x;
fo = fx;
if k == MaxIter
fprintf('Just best in %d iterations', MaxIter);
end
2.2 函数测试
文件名称:test.m
% 用最速下降法求最优化解
clc;
% 目标函数
fun = inline('x(1)^2+2*x(2)^2-2*x(1)*x(2)-4*x(1)', 'x');
% 目标函数的梯度函数
grad = inline('[2*x(1)-2*x(2)-4, 4*x(2)-2*x(1)]', 'x');
x0 = [0 0];
TolX = 1e-4;
TolFun = 1e-9;
MaxIter = 1000;
dist0 = 1;
[xo, fo] = Opt_Steepest(fun, grad, x0, TolX, TolFun, dist0, MaxIter)
2.3 测试结果
3 C程序
3.1 函数说明文件名称:function.c #include <math.h>
#include <stdio.h>
#include <stdlib.h>
# define eps 1e-6
double f (int coe[], double x[])
{
return (double) coe[0]*pow(x[0],2)+coe[1]*pow(x[1],2)+coe[2]*x[0]*x[1]+coe[3]*x[0]+coe[4]*x[1];
}
void grads (int coe[], double x[],double grads_x[])
{
grads_x[0] = (double) 2*coe[0]*x[0]+coe[2]*x[1]+coe[3];
grads_x[1] = (double) 2*coe[1]*x[1]+coe[2]*x[0]+coe[4];
}
double findLambda (int c[], double x[], double g[])
{
double inf_x[2] , lamd = 1;
do
{
inf_x[0] = -2*c[0]*(x[0]-lamd*g[0])*g[0]-2*c[1]*(x[1]-lamd*g[1])*g[1]-x[0]*g[1]*c[2]-x[1]*g[0]*c[2]+2*c[2]* g[0]*g[1]*lamd-g[0]*c[3]-c[4]*g[1];
inf_x[1] = 2*c[0]*g[0]*g[0]+2*c[1]*g[1]*g[1]+2*c[2]*g[0]*g[1];
lamd = lamd - inf_x[0]/inf_x[1];
}while(inf_x[0]>eps);
return lamd;
}
3.2 函数测试
文件名称:runme.c
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include "function.c"
# define eps 1e-6
int main()
{
int coe[5] = {1, 2, -2, -4, 0};
double x[2] = {0, 0},grads_x[2],lambda;
double fx = 0;
int i=0;
grads(coe, x, grads_x); /* 求梯度*/
while(sqrt(pow(grads_x[0],2)+pow(grads_x[1],2))-eps>0)
{
lambda = findLambda(coe, x, grads_x); /* 一维精确搜索*/
x[0] = x[0] - lambda*grads_x[0];
x[1] = x[1] - lambda*grads_x[1];
grads(coe, x, grads_x);
}
fx = f(coe, x);
printf("The optimal solution is: \n");
for(i=0;i<2;i++)
{
printf("%f\n",x[i]);
}
printf("\nThe optimal value is: ");
printf("\nfx = %f \n",fx);
getch();
return 0;
}
3.3 测试结果。