运筹学第六次实验报告

运筹学第六次实验报告

题目：用matlab编程最速下降法SteepestDescent 用matlab编程Newton法Newton

一、Armijo算法的算法框架：

1.给定初始点x0,x0?R，且令k=0；给定方向d; 2.给定初始步长?，通常取??1；

3.判断退出条件是否成立，若是则退出，否则执行step4，

退出条件：f(xk?1)?f(xk)?????[?f(xk)]T?dk?1,??(0.5,1)，通常??0.9； 4.

n?k?1????k,?为压缩系数;xk?1?xk??k?1?dk?1，返回step3 。

n二、最速下降法的算法框架：

1.给定初始点x0,x0?R，令k=0；

2.判断退出条件是否满足，若是则退出，否则执行step3，

退出条件：?f(xk)??；

3.选取方向dk?1，满足?f(xk)T?dk?1?0，最速下降法中dk?1???f(xk)； 4.利用Armijo算法确定步长5.令xk?1?xk??k?1；

?k?1?dk?1，k=k+1，返回step2。

n三、Newton法的算法框架：

1.给定初始点x0,x0?R，令k=0；

2.判断退出条件是否满足，若是则退出，否则执行step3，

退出条件：?f(xk)??；

3.选取方向dk?1，满足?f(xk)T?dk?1?0，Newton法中dk?1??H?1?f(xk), 其中H??f(xk)； 4.利用Armijo算法确定步长5.令xk?1?xk?2?k?1；

?k?1?dk?1，k=k+1，返回step2。

四、程序代码：

% Armijo方法获取非精确步长 function alpha = Armijo(x,fx,d,g,nf) % 输入--> x:迭代点 % 输入 x----迭代点

% fx---测试函数在x处的函数值 % d----下降方向

% g----测试函数在x处的梯度值 % nf---第nf个测试函数 % 输出 alpha----步长 mu = 0.9;

beta = 0.5; % 步长压缩系数 alpha = 1; % 初始步长 zeta = mu*g'*d;

temp = fx + alpha * zeta; xnew = x + alpha*d; fnew = fun(xnew,nf); iter = 0;

while fnew >= temp alpha = alpha*beta;

temp = fx + alpha * zeta; xnew = x + alpha*d; fnew = fun(xnew,nf); end end

% 最速下降法

function [x,fmin,iter]=SteepestDescent(x0,nf) % 输入 x0 ---- 初始点

% nf ---- 第nf个测试函数 % 输出 x ----- 极小值点

% fmin ---- 极小值，即x处的函数值 % iter ---- 迭代次数 epsilon = 1e-5; itermax = 1e4; [f,g] = fun(x0,nf); if isnan(g)

error('该函数没有梯度或未给出梯度！'); end

iter = 0; % 记录迭代次数

while norm(g)>epsilon && iter

alpha = Armijo(x0,f,d,g,nf); x0 = x0 + alpha * d; [f,g] = fun(x0,nf); iter = iter + 1; end x = x0; fmin = f; end

% Newton法

function [x,fmin,iter]=Newton(x0,nf) % 输入 x0 ---- 初始点

% nf ---- 第nf个测试函数 % 输出 x ----- 极小值点

% fmin ---- 极小值，即x处的函数值 % iter ---- 函数值计算次数 epsilon = 1e-5; itermax = 1e4;

[f,g,H] = fun(x0,nf); if isnan(g)

error('该函数没有梯度！'); end

if isnan(H)

error('该函数没有二阶梯度,不能用Newton法计算！'); end

iter = 0; % 记录迭代次数

while norm(g)>epsilon && iter

alpha = Armijo(x0,f,d,g,nf); x0 = x0 + alpha * d; [f,g,H] = fun(x0,nf); iter = iter + 1; end x = x0; fmin = f; end

function [f,g,H] = fun(x,nf)

% 计算指定第nf个函数在x点处的函数值f、一阶梯度值g和二阶梯度值H % 只计算二维例子 switch nf case 1

% p179 例6

% 最小值点(1,1), 最小值0

% 建议初始点(0 0)'

f = (x(1)-1)^2+(x(2)-1)^2;

g = zeros(2,1);

g(1,1) = 2*(x(1)-1); g(2,1) = 2*(x(2)-1);

H = zeros(2,2);

H(1,1) = 2; H(1,2) = 0; H(2,1) = 0; H(2,2) = 2; case 2

% Rosenbrock function

% 最小值点(1,1), 最小值0 % 建议初始点(-1 0.5)'

f = 100*(x(1)^2-x(2))^2+(x(1)-1)^2;

g = zeros(2,1);

g(1,1) = 400*x(1)*(x(1)^2-x(2))+2*(x(1)-1); g(2,1) = -200*(x(1)^2-x(2));

H = zeros(2,2);

H(1,1) = 400*(3*x(1)^2-x(2))+2; H(1,2) = -400*x(1); H(2,1) = H(1,2); H(2,2) = 200; otherwise

% modified Broyden tridiagonal function

% 来自经典文献Kolda-Lewis-Torczon-2003的389页例子 % 该例子不存在导数

% 最小值点(0,0), 最小值1 % 建议初始点(-0.9 -1)'

f = abs((3-2*x(1))*x(1)-2*x(2)+1)^(7/3) + abs((3-2*x(2))*x(2)-x(1)+1)^(7/3);

g = NaN; H = NaN; end

五、运行结果：

>> [x,fmin,iter]=SteepestDescent([-2,2]',2) x =

1.0000 1.0000

fmin =

1.1935e-010

iter =

2202

>> [x,fmin,iter]=Newton([-2,2]',2) x =

1.0000 1.0000

fmin =

1.6258e-012

iter =

196

由上述可知，newton法要比最速下降法要快的多。