无约束最优化问题之Rosenbrock法

算法目的

  求解无约束优化问题的极小值（最小值）。

算法特点

  算法步骤中不需要计算目标函数的导数。

算法步骤

例题1

  用Rosenbrock方法解
$$\min f(x):=(x_1-3{)}^2+2(x_2+2{)}^2.$$

解：设置初始基点坐标$x^{(0)}=(0,0{)}^T$，初始搜索方向集 d = { ( 1 , 0 ) T ,   ( 0 , 1 ) T } d=\{(1,0)^T,\ (0,1)^T \} d={ (1,0)T, (0,1)T}，初始步长集 δ = { 1 ,   1 } \delta =\{1,\ 1\} δ={ 1, 1}，加速因子$\alpha =3$，衰减因子$\beta =-0.5$，允许误差$ϵ=0.01$。

import numpy as np
import sys


# 定义目标函数
def function(x):
    return (x[0] - 3.) ** 2 + 2 * (x[1] + 2) ** 2


# 输入算法参数: 初始点坐标xk, 初始搜索方向sd, 初始步长delta, 加速因子alpha, 衰减因子beta, 允许误差epsilon
xk = np.array([0., 0])
sd = np.array([[1., 0], [0, 1]])
delta_0 = np.array([1., 1])
alpha = 3.
beta = -0.5
epsilon = 0.01

# 初始化大探测阶段数, 计算变量维数, 初始化第1大阶第1轮探测出发点和生成迭代用步长
k1, k2 = 1, 1
dim = len(xk)
y1 = xk.copy()
delta = delta_0.copy()

while 1:

    # 进入一个新的大探测阶段
    print('\n进入第', k1, '个大探测阶段')
    print('基点为:', xk, '函数值:', function(xk))
    print('各搜索方向:\n', sd)

    # 初始化用于构造新正交方向列的Lambda列
    Lambda = np.zeros([dim, 1])

    while 1:
        # 进入每一轮小探测
        print('\n  进入第', k2, '轮小探测')
        print('  探测出发点为:', y1)
        print('  各搜索步长:\n', delta)
        yk = y1.copy()
        for i in range(dim):
            t1, t2 = function(yk + delta[i] * sd[i]), function(yk)
            if t1 < t2:
                Lambda[i, 0] += delta[i]
                yk += delta[i] * sd[i]
                delta[i] *= alpha
            else:
                delta[i] *= beta
            print('  第', i + 1, '次探测得到的点;', yk, '函数值', function(yk))
        k2 += 1
        if function(yk) < function(y1):
            y1 = yk
        else:
            if function(yk) < function(xk):
                break
            elif False not in (np.abs(delta) <= epsilon):
                print('\n满足delta<=epsilon, 停止迭代, 最优解的估计:', xk, '函数值', function(xk))
                sys.exit()
            else:
                y1 = yk

    # 进入一个新的方向构造阶段
    if np.linalg.norm(yk - xk, ord=2) <= epsilon:
        print('\n满足||x(k+1)-x(k)|| <= epsilon, 停止迭代, 最优解的估计', yk, '函数值', function(yk))
        sys.exit()
    else:
        xk = yk
        y1 = yk

        # 新定义一组方向{p(n)}
        for i in range(dim):
            if Lambda[i] != 0:
                sd[i, :] = np.sum(Lambda[i:, 0] * sd[i:, :], axis=0)

        # 使用Gram-Schmidt正交化方法把{p(n)}正交化
        for i in range(1, dim):
            t = 0
            for j in range(i):
                t += (np.dot(sd[j, :], sd[i, :].T) / np.dot(sd[j, :], sd[j, :].T)) * sd[j, :]
            sd[i, :] -= t

        # 再单位化向量
        for i in range(dim):
            sd[i, :] /= np.linalg.norm(sd[i, :], ord=2)

        # 初始化步长
        delta = delta_0.copy()
    k1 += 1
    k2 = 1

  结果：

进入第 1 个大探测阶段
基点为: [0. 0.] 函数值: 17.0
各搜索方向:
 [[1. 0.]
 [0. 1.]]
 
  进入第 1 轮小探测
  探测出发点为: [0. 0.]
  各搜索步长:
 [1. 1.]
  第 1 次探测得到的点; [1. 0.] 函数值 12.0
  第 2 次探测得到的点; [1. 0.] 函数值 12.0
  
  进入第 2 轮小探测
  探测出发点为: [1. 0.]
  各搜索步长:
 [ 3