Python和Java代码实现：进退法确定搜索区间

引言

之前的文章中介绍过，黄金分割法可以求解如下的一维最优化问题：在一个搜索区间内，包含了目标函数的极小值点。如果此时只有函数的一个初始点$x_0$，但是并不知道极小值所在搜索区间，那么可以先使用本文即将介绍的进退法来确定搜索区间，然后再使用黄金分割法得到极值。

进退法原理

进退法的理论依据是：假定$f(x)$为单峰函数，只有一个极小值，且区间$[a,b]$为极小值的一个搜索区间，那么对于任意$x_1,x_2\in [a,b]$，如果$f(x_1)<f(x_2)$，则搜索区间可以收缩为$[a,x_2]$；如果$f(x_1)>f(x_2)$，则搜索区间可以收缩为$[x_1,b]$。

因此，如果给定初始点$x_0$和搜索步长$h$的情况下，首先以初始步长向前搜索一步，并计算对应的值$f(x_0+h)$:

（1）如果$f(x_0)<f(x_0+h)$，可以反推出搜索区间为$[\tilde{x},x_0+h]$，其中$\tilde{x}$虽然暂时未知，但是右端点可以确定为$x_0+h$；此时再后退一步计算$f(x_0-\lambda h)$，直到找到一个$\lambda$，使得$f(x_0-\lambda h)>f(x_0)$，从而确定搜索区间为$[f(x_0-\lambda h),f(x_0+h)]$。

（2）如果$f(x_0)>f(x_0+h)$，可以反推出搜索区间为$[x_0,\tilde{x}]$，其中$\tilde{x}$待求，此时前进一步计算$f(x_0+\lambda h)$，直到找到一个$\lambda$，使得$f(x_0+\lambda h)>f(x_0+h)$，从而确定搜索区间为$[f(x_0),f(x_0+\lambda h)]$。

（3）如果$f(x_0)=f(x_0+h)$，因为$f(x)$是单峰函数，所以搜索区间可以直接确定为$[f(x_0),f(x_0+h)]$。

python代码实现

一维最优化问题还是用黄金分割法中的实例：
$$min(t^2-5t+8)$$
该问题最小值为1.75，其对应的$t=2.5$。

实现以上进退法的Python代码如下所示。

# 待优化函数
def f(t):
    return t ** 2 - t * 5 + 8


# 进退法：确定搜索区间
def advance_and_retreat(func, x, h):

    if abs(func(x) - func(x + h)) <= 1e-6:
        # 第三种情况
        x_min, x_max = x, x + h
    elif func(x) < func(x + h):
        # 第一种情况
        x_max = x + h
        lamb = 1
        while func(x - lamb * h) < func(x):
            lamb += 1
        x_min = x - lamb * h
    else:
        # 第二组情况
        x_min = x + h
        lamb = 2
        while func(x + lamb * h) < func(x + h):
            lamb += 1
        x_max = x + lamb * h

    return x_min, x_max


if __name__ == '__main__':

    # 步长
    step = 0.1

    # 第一种情况实例
    x1 = 2
    lb1, ub1 = advance_and_retreat(f, x1, step)
    print('S1: x0 = {}, lb = {}, ub = {}'.format(x1, lb1, ub1))

    # 第二种情况实例
    x2 = 4
    lb2, ub2 = advance_and_retreat(f, x2, step)
    print('S2: x0 = {}, lb = {}, ub = {}'.format(x2, lb2, ub2))

    # 第三种情况实例
    x3 = 2.45
    lb3, ub3 = advance_and_retreat(f, x3, step)
    print('S3: x0 = {}, lb = {}, ub = {}'.format(x3, lb3, ub3))

主函数中，选取了3种不同的初值，分别对应了原理中所描述的3种情况。运行以上代码后，可以得到以下结果。显然每一个区间内都包含了极值点2.5，即都得到了正确的搜索区间。

S1: x0 = 2, lb = 2.1, ub = 2.9
S2: x0 = 4, lb = 1.0, ub = 4.1
S3: x0 = 2.45, lb = 2.45, ub = 2.5500000000000003

Java代码实现

Java代码如下所示。

public class AdvanceAndRetreat {
    public static void main(String[] args) {
        // 步长
        double step = 0.1;

        // 第一种情况实例
        double x1 = 2.0;
        Intervals interval1 = advanceAndRetreatMethod(x1, step);
        System.out.println("S1: x0 = " + x1 + ", lb = " + interval1.lb +  ", ub = " + interval1.ub);

        // 第二种情况实例
        double x2 = 4;
        Intervals interval2 = advanceAndRetreatMethod(x2, step);
        System.out.println("S2: x0 = " + x2 + ", lb = " + interval2.lb +  ", ub = " + interval2.ub);

        // 第三种情况实例
        double x3 = 2.45;
        Intervals interval3 = advanceAndRetreatMethod(x3, step);
        System.out.println("S2: x0 = " + x3 + ", lb = " + interval3.lb +  ", ub = " + interval3.ub);

    }

    // 进退法：确定搜索区间
    private static Intervals advanceAndRetreatMethod(double x, double h) {

        double lb = 0, ub = 0;

        if (Math.abs(func(x) - func(x+h)) <= 1e-6) {
            // 第三种情况
            lb = x;
            ub = x + h;
        }
        else if (func(x) < func(x+h)) {
            // 第一种情况
            ub = x + h;
            int lamb = 1;
            while (func(x-lamb * h) < func(x)) {
                lamb += 1;
            }
            lb = x - lamb * h;
        }
        else {
            // 第二种情况
            lb = x + h;
            int lamb = 2;
            while (func(x+lamb * h) < func(x + h)) {
                lamb += 1;
            }
            ub = x + lamb * h;
        }

        // 构造搜索区间对象
        Intervals interval = new Intervals();
        interval.lb = lb;
        interval.ub = ub;

        return interval;

    }

    // 待优化函数
    private static double func(double t) {
        return t * t - t * 5 + 8;
    }

    // 搜索区间对象
    private static class Intervals {
        double lb;
        double ub;
    }

}

运行以上代码后，可以得到

S1: x0 = 2.0, lb = 2.1, ub = 2.9
S2: x0 = 4.0, lb = 1.0, ub = 4.1
S2: x0 = 2.45, lb = 2.45, ub = 2.5500000000000003