matlab实现进退法、0.618法、牛顿法

实验一 一维最优化方法

实验目的及要求：

1. 编写程序实现进退法，利用进退法求解极值区间实例；

2. 编写程序实现0.618 法，利用0.618 法求解极值实例；

3. 编写程序实现牛顿法，来求解一维函数的局部极小值点。

实验原理：

在$Matlab$ 环境下，按照要求编写函数和程序，求解实例，直至取得正确的运行结果。（写出每种算法的步骤）

算法1：进退法

1. 设定初始值$\alpha =2，k=0，{\varphi }_0=\varphi ({\lambda }_0)$，传入参数$h_0$。
2. 使${\lambda }_{k+1}={\lambda }_k+h_k$，计算${\varphi }_{k+1}=\varphi ({\lambda }_{k+1})$。
3. 若${\varphi }_{k+1}<{\varphi }_k$，转向4；反之转向5。
4. 加步，$h_{k+1}=a{h}_k，\lambda ={\lambda }_k，{\lambda }_k={\lambda }_{k+1}，{\varphi }_k={\varphi }_{k+1}，k++$，转向2。
5. 反向探索：
   • 若$k=0$则转换方向：$h_k=-h_k，\lambda ={\lambda }_{k+1}$，转向2。
   • 否则停止迭代，$a=min(\lambda ,{\lambda }_{k+1})，b=max(\lambda ,{\lambda }_{k+1})$，输出$[a,b]$。

算法2：0.618法

1. 选取初始区间$[a_1,b_1]$和误差$\epsilon >0$，$\tau$=0.618。
2. 计算${\lambda }_1=a_1+(1-\tau )\left(b_1-a_1\right),{\mu }_1=a_1+\tau \left(b_1-a_1\right)$，并求出$\varphi ({\lambda }_1)，\varphi ({\mu }_1)$，$k=1$。
3. 若$\left|b_k-a_k\right|<\epsilon$，停止计算，输出$\bar{t}=\frac{{\lambda }_k+{\mu }_k}{2}$；反之转向4。
4. 比较函数值
   • $\varphi ({\lambda }_k)<\varphi ({\mu }_k)$，转向5。
   • 反之，转向6。
5. 向左搜索，令$a_{k+1}=a_k,b_{k+1}={\mu }_k,{\mu }_{k+1}={\lambda }_k,\phi \left({\mu }_{k+1}\right)=\phi \left({\lambda }_k\right)$，并计算${\lambda }_{k+1}=a_{k+1}+(1-\tau )\left(b_{k+1}-a_{k+1}\right)$和$\phi \left({\lambda }_{k+1}\right)$，转向7。
6. 向右搜索，令$a_{k+1}={\lambda }_k,b_{k+1}=b_k,{\lambda }_{k+1}={\mu }_k,\phi \left({\lambda }_{k+1}\right)=\phi \left({\mu }_k\right)$，并计算${\mu }_{k+1}=a_{k+1}+\tau \left(b_{k+1}-a_{k+1}\right)$和$\phi \left({\mu }_{k+1}\right)$，转向7。
7. $k++$，转向3。

算法3：牛顿法

1. 设定初始点$t_0$，允许误差$\epsilon >0$，$k=0$。
2. ${\varphi }^{\prime }(t_k)<\epsilon$输出$t_k$，停止计算，反之转向3。
3. 计算$t_{k+1}=t_k-\frac{{\phi }^{\prime }\left(t_k\right)}{{\phi }^{\prime \prime }\left(t_k\right)}，k++$，转向2。

实验内容（方法和步骤）

例题

题目1 编写程序实现进退法。

利用Matlab 编写函数 [xmin,xmax]=Brackeing(f,x0,h0). 区间[xmin,xmax] 包含函数f的一个局部极值点。

调用上述函数，实现代码在代码实现部分：

function [xmin,xmax]=Brackeing(f,x0,h)
%{
Brackeing 进退法
   f       函数
   x0      初始点
   h0      步长
   
   xmin    左区间
   xmax    右区间
%}
F = matlabFunction(f);  %建立函数句柄
a = 2;                  %加倍系数
k = 1;                  %迭代数
lambda = 0;             %记录答案
flag = false;           %标记

x = zeros(1,10);        %申请自变量空间
f_val = zeros(1,10);    %申请函数值空间
x(1) = x0;              %自变量赋值
f_val(1) = F(x(1));     %函数赋值


while true
%%%%%%%%%%%%% 试探&计算
    x(k+1) = x(k) + h;      
    f_val(k+1) = F(x(k+1));
%%%%%%%%%%%%%
    if f_val(k+1) >= f_val(k) %  比较
        if k == 1       %检查
            h = -h;     %改变方向反向前进
            lambda = x(k+1);
            x(k+1) = x(k) + h;      
            f_val(k+1) = F(x(k+1));
        else 
            flag = true;
        end
    end
    
    if flag         %结束计算
        break
    end
        
%%%%%%%%%%%%% 加步
    h = a * h;
    lambda = x(k);
    x(k) = x(k+1);
    f_val(k) = f_val(k+1);
    k = k + 1;
%%%%%%%%%%%%%
end

xmin = min(lambda, x(k+1));     %左边界
xmax = max(lambda, x(k+1));     %右边界
disp(xmin)
disp(xmax)
%%%%%%%%%%%%% 绘图查看
figure('Name', 'Brackeing')
xx = xmin-0.5:0.01:xmax+0.5;
yy = F(xx);
plot(xx, yy, 'LineWidth', 2)
hold on
plot([xmin,xmax], [min(yy), min(yy)],':r', 'LineWidth', 1.5)
end

题目2 利用进退法求解极值区间实例

利用进退法求解极值区间实例。取初始点x0=0，步长h0=0.1 ，用进退法求函数$f(x)=(x^2-1{)}^2-(x-1{)}^2$+3的极值区间。

代码实现部分为：

syms x y;
y = (x^2 - 1)^2 - (x - 1)^2 + 3;
[xmin,xmax]=Brackeing(y,0,0.1);

输出结果为：

xmin = -3.1000
xmax = -0.7000

检查结果为：

题目3 编写程序实现0.618 法

利用Matlab 编写函数 . 调用函数得到的x 为函数f 在区间[a,b] 上的局部极小值点。

代码实现部分为：

function [x,miny]=MinHJ(f,a,b,eps)
%{
    MinHJ   0.618法
    f       函数
    a       左区间
    b       右区间
    eps     精度

    x       自变量
    miny    局部极小值
%}
xmax = b;
xmin = a;
F = matlabFunction(f);  %函数句柄
tau = 0.618;            %缩短率
lambda = a + (1 - tau)*(b - a);     %计算lambda
mu = a + tau * (b - a);             %计算mu

lambda_val = F(lambda);             %lambda函数值
mu_val = F(mu);                     %mu函数值

while b - a >= eps                  %误差区间之外
    %fprintf('a=%d,b=%d,a_val=%d,b_val=%d\n',a,b,lambda_val,mu_val);
    if lambda_val < mu_val
        b = mu;
        mu = lambda;
        mu_val = lambda_val;
        lambda = a + (1 - tau)*(b - a);
        lambda_val = F(lambda);
    else
        a = lambda;
        lambda = mu;
        lambda_val = mu_val;
        mu = a + tau * (b - a);
        mu_val = F(mu);
    end
end

%计算结果
x = (lambda + mu) / 2;
miny = F(x);

%%%%%%%% 绘图
figure('Name', 'MinHJ')
xx = xmin:0.01:xmax;
yy = F(xx);
plot(xx, yy, 'LineWidth', 2)
hold on
plot(x,miny,'p', 'LineWidth', 2)
end

题目4 利用0.618 法求解极值实例

利用0.618 法求下面函数的极小值:$f(x)={\left(x^2-1\right)}^2-(x-1{)}^2+3,t\in (-10,10)$，误差为0.01

调用上述代码，实现部分为：

syms x y;
y = (x^2 - 1)^2 - (x - 1)^2 + 3;
[x,miny]=MinHJ(y,-10,10,0.01);

输出结果为：

x = -1.3656
miny = -1.8481

检查结果为：

题目5 编写程序实现牛顿法

利用Matlab 编写函数 [x,minf]=Newton(f,x0,eps). 其中f 为目标函数， x0 为初始点，eps 为算法终止的精度。用牛顿法求解无约束优化问题minf(x) 。

代码实现部分为：

function [x,minf]=Newton(f,x0,eps)
%{
Newton      牛顿法
    f       函数
    x0      初始点
    eps     误差
   
    x       目标x
    minf    目标函数值
%}
xmax = x0+2;
xmin = x0-2;
F0 = matlabFunction(f);             %原函数函数句柄
F1 = matlabFunction(diff(f));       %一阶导函数函数句柄
F2 = matlabFunction(diff(f,2));     %二阶导函数函数句柄

while abs(F1(x0)) >= eps            %迭代
    x0 = x0 - F1(x0) / F2(x0);
end

x = x0;
minf = F0(x);

%%%%%%%% 绘图
figure('Name', 'MinHJ')
xx = xmin:0.01:xmax;
yy = F0(xx);
plot(xx, yy, 'LineWidth', 2)
hold on
plot(x,minf,'p', 'LineWidth', 2)
end


while abs(F1(x0)) >= eps            %迭代
    x0 = x0 - F1(x0) / F2(x0);
end

x = x0;
minf = F0(x);

%%%%%%%% 绘图
figure('Name', 'MinHJ')
xx = xmin:0.01:xmax;
yy = F0(xx);
plot(xx, yy, 'LineWidth', 2)
hold on
plot(x,minf,'p', 'LineWidth', 2)
end