[最优化技术] 3-7 坐标轮换法

3-7 坐标轮换法

坐标轮换法

简述

\(\qquad\)上一篇介绍了修正牛顿法；本篇转向一种不需要梯度信息的直接搜索方法。坐标轮换法（Coordinate Descent Method），又称变量轮换法、交替法或降维法。

\(\qquad\)与最速下降法类似，它也是一种迭代搜索方法。不同之处在于，最速下降法使用的是负梯度方向，而坐标轮换法直接使用坐标轴方向（基向量）作为搜索方向。

\(\qquad\)该方法的基本思想是：每次只沿一个坐标轴方向进行一维搜索，找到该方向上的极小值点后，再切换到下一个坐标轴方向进行搜索。轮流沿各个坐标轴方向搜索，直到满足收敛条件。这种方法将多维问题转化为一系列的一维问题，因此也被称为“降维法”。

原理

对于 \(n\) 元函数 \(f(X)\)，其搜索方向取为坐标系的基向量 \(e_i\)：

\[S_i^{(k)} = e_i = (0, \dots, 1, \dots, 0)^T \quad (1 \text{ 在第 } i \text{ 个位置}) \]

在一轮迭代中，算法依次沿 \(e_1, e_2, \dots, e_n\) 方向进行 \(n\) 次一维搜索。第 \(k\) 轮第 \(i\) 次搜索的迭代式为

\[X_i^{(k)} = X_{i-1}^{(k)} + \alpha_i^{(k)} S_i^{(k)}, \quad i=1,2,\dots,n \]

其中 \(\alpha_i^{(k)}\) 为沿 \(S_i^{(k)}\) 方向作一维搜索得到的最优步长。

注意：坐标轮换法的搜索方向是固定的，\(S_i^{(k)}\) 与迭代轮数 \(k\) 无关。

具体步骤

对于给定初始点 \(X_0^{(0)}\) 和收敛精度 \(\varepsilon\)：

① 确定初始点 \(X_0^{(0)}\)，收敛精度 \(\varepsilon\)，令轮数 \(k=1\)；

② 取基向量作为搜索方向，\(S_i^{(k)} = e_i, \; i=1, 2, \dots, n\)（对于二元函数 \(n=2\)）；

③ 从当前点 \(X_0^{(k)}\) 出发：

  • 沿 \(S_1^{(k)}\) 方向进行一维搜索，得到点 \(X_1^{(k)}\)；
  • 从 \(X_1^{(k)}\) 出发，沿 \(S_2^{(k)}\) 方向进行一维搜索，得到点 \(X_2^{(k)}\)；
  • \(\dots\)
  • 从 \(X_{n-1}^{(k)}\) 出发，沿 \(S_n^{(k)}\) 方向进行一维搜索，得到本轮终点 \(X_n^{(k)}\)。

即：\(X_0^{(k)} \xrightarrow{S_1} X_1^{(k)} \xrightarrow{S_2} X_2^{(k)} \dots \xrightarrow{S_n} X_n^{(k)}\)。

④ 收敛判断：若 \(\|X_n^{(k)} - X_0^{(k)}\| \le \varepsilon\)，则结束，最优解为 \(X_n^{(k)}\)；否则，令 \(X_0^{(k+1)} = X_n^{(k)}\)，\(k=k+1\)，重复步骤③。

符号说明：

  • \(X_i^{(k)}\) 表示：在第 \(k\) 轮迭代中，沿第 \(i\) 个坐标轴方向进行一维搜索后得到的点。它既是第 \(i\) 次搜索的终点，也是第 \(i+1\) 次搜索的起点。
  • 对于 \(n\) 元函数，每 1 轮迭代都需要进行 \(n\) 次一维搜索。

程序框图

坐标轮换法的程序框图

符号说明：图中 \(x\)、\(d_i^k\)、\(\alpha_i^k\) 分别对应正文的 \(X\)、\(S_i^{(k)}\)、\(\alpha\)；图中轮数从 \(k=0\) 起算，判据为 \(<\varepsilon\)；正文从 \(k=1\) 起算，判据为 \(\le\varepsilon\)。

实例计算

这里沿用之前的例子，演示前两轮迭代，并给出收敛过程。

题目

\[f(X) = x_1^2 + 2x_2^2 - 4x_1 - 2x_1 x_2,\;\; X^{(0)} = [1,1]^T,\;\varepsilon = 0.01 \]

记 \(X_0^{(1)} = X^{(0)}\)。

解答

为书写方便，下面把每次一维搜索中被优化的那个变量统一记为 \(\alpha\)（它在该次搜索中就是自由坐标，不是上文的最优步长），对应的一元函数记为 \(\varphi(\alpha)\)。

第 1 轮迭代 (\(k=1\))：

初始点 \(X_0^{(1)} = [1, 1]^T\)。

1) 沿 \(x_1\) 轴方向搜索 (\(S_1 = [1, 0]^T\))：

固定 \(x_2 = 1\)，函数变为关于 \(\alpha\) 的一元函数：

\[\varphi(\alpha) = f(\alpha, 1) = \alpha^2 + 2(1)^2 - 4\alpha - 2\alpha(1) = \alpha^2 - 6\alpha + 2 \]

求导令其为 0：\(\varphi'(\alpha) = 2\alpha - 6 = 0 \implies \alpha = 3\)。
得到新点 \(X_1^{(1)} = [3, 1]^T\)。

2) 沿 \(x_2\) 轴方向搜索 (\(S_2 = [0, 1]^T\))：

从 \(X_1^{(1)} = [3, 1]^T\) 出发，固定 \(x_1 = 3\)，同样用 \(\alpha\) 表示当前自由变量：

\[\varphi(\alpha) = f(3, \alpha) = 3^2 + 2\alpha^2 - 4(3) - 2(3)\alpha = 2\alpha^2 - 6\alpha - 3 \]

求导令其为 0：\(\varphi'(\alpha) = 4\alpha - 6 = 0 \implies \alpha = 1.5\)。
得到本轮终点 \(X_2^{(1)} = [3, 1.5]^T\)。

3) 收敛判断：

\[\|X_2^{(1)} - X_0^{(1)}\| = \|[3, 1.5]^T - [1, 1]^T\| = \|[2, 0.5]^T\| = \sqrt{4 + 0.25} \approx 2.06 > \varepsilon \]

不满足收敛条件，进入下一轮。

第 2 轮迭代 (\(k=2\))：

从 \(X_0^{(2)} = X_2^{(1)} = [3, 1.5]^T\) 出发。

1) 沿 \(x_1\) 轴方向搜索 (\(S_1 = [1, 0]^T\))：

固定 \(x_2 = 1.5\)：

\[\varphi(\alpha) = f(\alpha, 1.5) = \alpha^2 + 2(1.5)^2 - 4\alpha - 2\alpha(1.5) = \alpha^2 - 7\alpha + 4.5 \]

求导得 \(\varphi'(\alpha) = 2\alpha - 7 = 0 \implies \alpha = 3.5\)，所以 \(X_1^{(2)} = [3.5, 1.5]^T\)。

2) 沿 \(x_2\) 轴方向搜索 (\(S_2 = [0, 1]^T\))：

固定 \(x_1 = 3.5\)：

\[\varphi(\alpha) = f(3.5, \alpha) = 2\alpha^2 - 7\alpha - 1.75 \]

求导得 \(\varphi'(\alpha) = 4\alpha - 7 = 0 \implies \alpha = 1.75\)，所以 \(X_2^{(2)} = [3.5, 1.75]^T\)。

3) 收敛判断：

\[\|X_2^{(2)} - X_0^{(2)}\| = \|[0.5, 0.25]^T\| = \sqrt{0.3125} \approx 0.56 > \varepsilon \]

仍不满足，继续迭代。

精确一维搜索下，各轮结果如下（第 1、2 轮的计算见上）：

轮数 \(k\) 本轮终点 \(X_2^{(k)}\) 位移 \(\lVert X_2^{(k)} - X_0^{(k)} \rVert\) \(f(X_2^{(k)})\)
1 (3.000000, 1.500000) 2.061553 -7.500000
2 (3.500000, 1.750000) 0.559017 -7.875000
3 (3.750000, 1.875000) 0.279508 -7.968750
4 (3.875000, 1.937500) 0.139754 -7.992188
5 (3.937500, 1.968750) 0.069877 -7.998047
6 (3.968750, 1.984375) 0.034939 -7.999512
7 (3.984375, 1.992188) 0.017469 -7.999878
8 (3.992188, 1.996094) 0.008735 -7.999969

第 8 轮位移 \(0.008735 \le \varepsilon\)，停止迭代，得到

\[X^* \approx [3.992188,\ 1.996094]^T,\qquad f(X^*) \approx -7.999969 \]

已经很接近理论最优解 \([4,2]^T\) 和 \(f^*=-8\)。

方法分析

  1. 方向简单，无需用梯度：外层搜索只需要函数值，不需要计算偏导数；一维搜索可用黄金分割法等只依赖函数值的方法，对导数难以求出的函数非常友好。

  2. 多维问题效率较低，收敛慢：每次只改变一个变量，其他变量固定，搜索路径呈直角锯齿状。等值线为长短轴不平行于坐标轴的椭圆时，需要多轮迭代，效率显著降低，如图 4-21(b) 所示。

  3. 特殊情况一轮到位：对于等高线为同心圆或长短轴与坐标轴平行的椭圆，即形如 \(F(X) = ax_1^2 + bx_2^2 + \text{线性项}\) 的函数（变量之间无交叉项 \(x_1 x_2\)），只需 1 轮迭代（\(n\) 次搜索）即可到达极值点；对二元函数就是两次一维搜索，如图 4-21(a) 所示。注意本节例题含 \(x_1 x_2\) 交叉项，等值线是斜置椭圆，属于 (b) 的情况，所以需要多轮迭代。

  4. 失效情况（脊线问题）：当等高线出现与坐标轴斜交的“脊线”时，沿坐标轴方向的搜索会停在脊线上，无法沿脊线方向继续前进，因而找不到最优点；这时需要进行转轴变换（坐标变换），如图 4-21(c) 所示。

坐标轮换法的效能

代码示例

下面只保留核心循环；函数定义、一维搜索子程序和输出见下方的完整程序。

C++代码示例（坐标轮换法核心循环）： 
// 核心循环：每轮依次沿 n 个坐标轴方向作一维搜索
for (int k = 0; k < max_rounds; ++k) {
    vector<double> x_start = x;

    for (int i = 0; i < n; ++i) {
        vector<double> d(n, 0.0);
        d[i] = 1.0;

        auto phi = [&](double alpha) {
            vector<double> new_x(x.size());
            for (size_t j = 0; j < x.size(); ++j)
                new_x[j] = x[j] + alpha * d[j];
            return func(new_x);
        };

        // 省略异常处理，完整写法见下
        auto bracket = bracket_minimum(phi, 0.0, 0.1);
        double alpha_star = golden_section(phi, bracket.first, bracket.second, eps);

        for (size_t j = 0; j < x.size(); ++j)
            x[j] += alpha_star * d[j];
    }

    double dist = 0.0;
    for (int i = 0; i < n; ++i)
        dist += (x[i] - x_start[i]) * (x[i] - x_start[i]);
    if (sqrt(dist) <= eps) break;
}

完整程序

C++代码实现（完整）： 
#include <iostream>
#include <vector>
#include <cmath>
#include <functional>
#include <iomanip>
#include <utility>

using namespace std;

// 进退法求初始区间 (复用)
pair<double, double> bracket_minimum(function<double(double)> func, double a0, double h) {
    const int MAX_ITER = 1e4;
    if (h <= 0) throw invalid_argument("Step size h must be positive.");
    if (func(a0) <= func(a0 + h)) {
        double x0 = a0 + h, x1 = a0;
        h *= 2;
        double x2 = x1 - h;
        int iter = 0;
        while (func(x2) < func(x1)) {
            if (iter++ > MAX_ITER) throw runtime_error("Max iterations exceeded");
            x0 = x1; x1 = x2; h *= 2; x2 = x1 - h;
        }
        return make_pair(x2, x0);
    } else {
        double x0 = a0;
        double x1 = a0 + h;
        h *= 2;
        double x2 = x1 + h;
        int iter = 0;
        while (func(x2) < func(x1)) {
            if (iter++ > MAX_ITER) throw runtime_error("Max iterations exceeded");
            x0 = x1; x1 = x2; h *= 2; x2 = x1 + h;
        }
        return make_pair(x0, x2);
    }
}

// 黄金分割法 (复用)
double golden_section(function<double(double)> func, double l, double r, double eps) {
    double a = l, b = r, aa, bb;
    while (b - a > eps) {
        aa = b - 0.618 * (b - a);
        bb = a + 0.618 * (b - a);
        if (func(aa) < func(bb)) b = bb;
        else a = aa;
    }
    return (a + b) / 2;
}

// 坐标轮换法
vector<double> coordinate_descent(
    function<double(const vector<double>&)> func,
    const vector<double>& x0,
    double eps,
    int max_rounds = 1000
) {
    vector<double> x = x0;
    int n = x0.size();
    
    cout << "\n迭代过程:" << endl;
    cout << "轮数\t点 X\t\t\t函数值 f(X)" << endl;
    cout << "--------------------------------------------------------" << endl;

    for (int k = 0; k < max_rounds; ++k) {
        vector<double> x_start = x;

        for (int i = 0; i < n; ++i) {
            vector<double> d(n, 0.0);
            d[i] = 1.0;

            auto phi = [&](double alpha) {
                vector<double> new_x(x.size());
                for (size_t j = 0; j < x.size(); ++j) new_x[j] = x[j] + alpha * d[j];
                return func(new_x);
            };

            pair<double, double> bracket;
            try { bracket = bracket_minimum(phi, 0.0, 0.1); } 
            catch (...) { bracket = make_pair(0.0, 10.0); }

            double alpha_star = golden_section(phi, bracket.first, bracket.second, eps);

            for (size_t j = 0; j < x.size(); ++j) x[j] += alpha_star * d[j];
        }

        cout << k + 1 << "\t";
        for (size_t i = 0; i < x.size(); ++i) {
            printf("%.6lf", x[i]);
            cout << (i == x.size() - 1 ? "" : ", ");
        }
        printf("\t%.6lf\n", func(x));

        double dist = 0.0;
        for (int i = 0; i < n; ++i) dist += (x[i] - x_start[i]) * (x[i] - x_start[i]);
        if (sqrt(dist) <= eps) {
            cout << "收敛！" << endl;
            break;
        }
    }
    return x;
}

int main() {
    cout << "坐标轮换法求解多元函数极小值" << endl;
    cout << "目标函数: f(X) = x1^2 + 2*x2^2 - 4*x1 - 2*x1*x2" << endl;
    cout << "理论最优解: X* = [4.0, 2.0]^T, f(X*) = -8.0" << endl << endl;

    auto f = [](vector<double> x) {
        return x[0] * x[0] + 2 * x[1] * x[1] - 4 * x[0] - 2 * x[0] * x[1];
    };

    vector<double> x0 = {1.0, 1.0};
    double eps = 0.01;

    vector<double> result = coordinate_descent(f, x0, eps);

    cout << "\n最终结果:" << endl;
    cout << "最优解 X* = [";
    for (size_t i = 0; i < result.size(); ++i) cout << result[i] << (i == result.size() - 1 ? "" : ", ");
    cout << "]^T" << endl;
    cout << "目标函数值 f(X*) = " << f(result) << endl;

    return 0;
}

\(\qquad\)程序用黄金分割法（\(\varepsilon=0.01\)）作一维搜索，实测第 1 轮为 \((3.002196,1.501746)\)，与手算的 \((3,1.5)\) 略有差异；第 8 轮收敛到 \((3.995681,1.997798)\)，\(f=-7.999991\)，末轮位移约 \(0.0097 \le 0.01\)。差异来自一维搜索的精度。

总结

坐标轮换法是一种最直接的降维优化方法。

优点：

  • 原理极其简单，无需计算梯度。
  • 每次一维搜索可达到给定精度；若用解析法，理论上可精确求解。
  • 对于等值线为圆或长短轴平行于坐标轴的椭圆的问题，效率很高，两次搜索即可到达最优点。

缺点：

  • 收敛速度通常较慢，且受等值线形状影响很大。
  • 等值线为斜置椭圆时，需要多轮迭代，效率显著降低。
  • 出现脊线时，搜索会停在脊线上，找不到最优点，需要转轴变换。

\(\qquad\)坐标轮换的逐维搜索思想也为后续方法提供了基础：鲍威尔法正是利用坐标轮换的原理构造搜索方向。

\(\qquad\)下一篇将介绍共轭方向法。它从坐标轮换法和最速下降法出发，构造一组共轭方向，使二次函数最多 \(n\) 步收敛；再下一篇介绍共轭梯度法。

本学习笔记参考资料：

[1] 白清顺, 孙靖民, 梁迎春. 机械优化设计 第7版[M]. 北京: 机械工业出版社, 2024. ISBN: 978-7-111-75103-8 在线链接

[2] 武汉理工大学《最优化技术B》课程课件，授课教师：颜彬老师

posted @ 2026-10-21 12:09  H_Elden  阅读(43)  评论(0)    收藏  举报