每日算法 — 使用java实现模拟退火:温度调度与全局最优解搜索

模拟退火(Simulated Annealing)是一种源于固体退火过程的通用概率型优化算法。在材料科学中,高温使金属原子自由运动,缓慢降温则让原子逐渐有序排列形成低能量晶格。模拟退火将这一物理过程抽象为数学模型:以一定概率接受劣解,跳出局部最优,随着”温度”降低逐步收敛到全局最优。本文将用Java完整实现模拟退火算法,并通过旅行商问题与函数优化两个经典案例展示其强大能力。

算法原理与Metropolis准则

模拟退火的核心在于Metropolis接受准则:对于当前解S和新解S’,计算能量差ΔE = E(S’) – E(S)。

  • 若ΔE < 0(新解更优),必定接受
  • 若ΔE ≥ 0(新解更劣),以概率exp(-ΔE/T)接受,其中T为当前温度

温度T随迭代逐渐降低,接受劣解的概率也随之减小。初期高温允许广泛探索,后期低温精细局部搜索,这种机制有效平衡了”探索”与”利用”。

温度调度策略

  1. 指数降温:T_{k+1} = α × T_k,α通常取0.95~0.995
  2. 对数降温:T_k = T_0 / ln(1+k),理论保证收敛但速度极慢
  3. 线性降温:T_{k+1} = T_k – ΔT,简单直接
  4. 自适应降温:根据接受率动态调整降温速度

Java实现

import java.util.*;

/**
 * 模拟退火算法通用框架
 * 支持任意解空间与能量函数,通过泛型实现高度可扩展
 */
public class SimulatedAnnealing<S> {

    /**
     * 解的邻居生成函数接口
     */
    @FunctionalInterface
    public interface NeighborGenerator<S> {
        S generate(S current);
    }

    /**
     * 能量函数接口(目标函数,越小越好)
     */
    @FunctionalInterface
    public interface EnergyFunction<S> {
        double calculate(S solution);
    }

    /**
     * 温度调度策略
     */
    @FunctionalInterface
    public interface CoolingSchedule {
        double nextTemperature(double currentTemp, int iteration);
    }

    private final NeighborGenerator<S> neighborGenerator;
    private final EnergyFunction<S> energyFunction;
    private final CoolingSchedule coolingSchedule;
    private final Random random;

    // 算法参数
    private double initialTemperature;
    private int maxIterations;
    private int iterationsPerTemp; // 每个温度下的迭代次数

    public SimulatedAnnealing(NeighborGenerator<S> neighborGenerator,
                               EnergyFunction<S> energyFunction,
                               CoolingSchedule coolingSchedule) {
        this.neighborGenerator = neighborGenerator;
        this.energyFunction = energyFunction;
        this.coolingSchedule = coolingSchedule;
        this.random = new Random();
    }

    /**
     * 设置算法参数
     */
    public void setParameters(double initialTemperature, int maxIterations, int iterationsPerTemp) {
        this.initialTemperature = initialTemperature;
        this.maxIterations = maxIterations;
        this.iterationsPerTemp = iterationsPerTemp;
    }

    /**
     * 执行模拟退火搜索
     * @param initialSolution 初始解
     * @return 搜索过程中的最优解
     */
    public SearchResult<S> search(S initialSolution) {
        S currentSolution = initialSolution;
        double currentEnergy = energyFunction.calculate(currentSolution);

        S bestSolution = currentSolution;
        double bestEnergy = currentEnergy;

        double temperature = initialTemperature;
        List<Double> energyHistory = new ArrayList<>();

        for (int iter = 0; iter < maxIterations; iter++) {
            for (int inner = 0; inner < iterationsPerTemp; inner++) {
                // 生成邻居解
                S neighbor = neighborGenerator.generate(currentSolution);
                double neighborEnergy = energyFunction.calculate(neighbor);
                double deltaE = neighborEnergy - currentEnergy;

                // Metropolis准则
                if (deltaE < 0 || random.nextDouble() < Math.exp(-deltaE / temperature)) {
                    currentSolution = neighbor;
                    currentEnergy = neighborEnergy;

                    if (currentEnergy < bestEnergy) {
                        bestEnergy = currentEnergy;
                        bestSolution = currentSolution;
                    }
                }
            }

            energyHistory.add(bestEnergy);

            // 降温
            temperature = coolingSchedule.nextTemperature(temperature, iter);

            // 提前终止条件:温度极低且长时间无改进
            if (temperature < 1e-10) {
                break;
            }
        }

        return new SearchResult<>(bestSolution, bestEnergy, energyHistory);
    }

    /**
     * 搜索结果封装
     */
    public static class SearchResult<S> {
        public final S bestSolution;
        public final double bestEnergy;
        public final List<Double> energyHistory;

        public SearchResult(S bestSolution, double bestEnergy, List<Double> energyHistory) {
            this.bestSolution = bestSolution;
            this.bestEnergy = bestEnergy;
            this.energyHistory = energyHistory;
        }
    }

    // 预定义降温策略
    public static CoolingSchedule exponentialCooling(double alpha) {
        return (temp, iter) -> temp * alpha;
    }

    public static CoolingSchedule linearCooling(double delta) {
        return (temp, iter) -> Math.max(0, temp - delta);
    }

    public static CoolingSchedule logarithmicCooling(double t0) {
        return (temp, iter) -> t0 / Math.log(2 + iter);
    }

    /**
     * 案例1:函数优化(Rastrigin函数)
     * Rastrigin是经典的多峰测试函数,大量局部极小值考验全局搜索能力
     */
    public static void functionOptimizationDemo() {
        System.out.println("=== 案例1: Rastrigin函数优化 ===");

        // Rastrigin函数:f(x) = 10n + Σ(xi² - 10cos(2πxi))
        final int dimensions = 10;
        final double searchRange = 5.12;

        EnergyFunction<double[]> rastrigin = (x) -> {
            double sum = 10 * dimensions;
            for (double xi : x) {
                sum += xi * xi - 10 * Math.cos(2 * Math.PI * xi);
            }
            return sum;
        };

        NeighborGenerator<double[]> neighborGen = (current) -> {
            double[] neighbor = current.clone();
            int idx = random.nextInt(dimensions);
            double step = (random.nextDouble() - 0.5) * 0.5; // 步长
            neighbor[idx] += step;
            // 边界约束
            neighbor[idx] = Math.max(-searchRange, Math.min(searchRange, neighbor[idx]));
            return neighbor;
        };

        // 初始解
        double[] initial = new double[dimensions];
        Random rand = new Random();
        for (int i = 0; i < dimensions; i++) {
            initial[i] = (rand.nextDouble() - 0.5) * 2 * searchRange;
        }

        double initialEnergy = rastrigin.calculate(initial);
        System.out.printf("初始能量: %.4f%n", initialEnergy);

        // 运行模拟退火
        SimulatedAnnealing<double[]> sa = new SimulatedAnnealing<>(
            neighborGen, rastrigin, exponentialCooling(0.99)
        );
        sa.setParameters(100.0, 5000, 10);

        SearchResult<double[]> result = sa.search(initial);
        System.out.printf("最优能量: %.4f%n", result.bestEnergy);
        System.out.printf("改进幅度: %.2f%%%n", (initialEnergy - result.bestEnergy) / initialEnergy * 100);
    }

    private static final Random random = new Random();

    /**
     * 案例2:旅行商问题(TSP)
     */
    public static void tspDemo() {
        System.out.println("\n=== 案例2: 旅行商问题(TSP) ===");

        // 生成随机城市坐标
        int cityCount = 30;
        double[][] cities = new double[cityCount][2];
        Random rand = new Random(42);
        for (int i = 0; i < cityCount; i++) {
            cities[i][0] = rand.nextDouble() * 100;
            cities[i][1] = rand.nextDouble() * 100;
        }

        // 计算两城市间距离
        double[][] distMatrix = new double[cityCount][cityCount];
        for (int i = 0; i < cityCount; i++) {
            for (int j = 0; j < cityCount; j++) {
                double dx = cities[i][0] - cities[j][0];
                double dy = cities[i][1] - cities[j][1];
                distMatrix[i][j] = Math.sqrt(dx * dx + dy * dy);
            }
        }

        // 能量函数:路径总长度
        EnergyFunction<int[]> tspEnergy = (path) -> {
            double total = 0;
            for (int i = 0; i < path.length - 1; i++) {
                total += distMatrix[path[i]][path[i + 1]];
            }
            total += distMatrix[path[path.length - 1]][path[0]]; // 回到起点
            return total;
        };

        // 邻居生成:2-opt交换(反转路径中一段)
        NeighborGenerator<int[]> tspNeighbor = (current) -> {
            int[] neighbor = current.clone();
            int i = random.nextInt(cityCount);
            int j = random.nextInt(cityCount);
            if (i > j) { int t = i; i = j; j = t; }
            // 反转i到j之间的路径段
            while (i < j) {
                int t = neighbor[i];
                neighbor[i] = neighbor[j];
                neighbor[j] = t;
                i++;
                j--;
            }
            return neighbor;
        };

        // 初始解:顺序访问
        int[] initialPath = new int[cityCount];
        for (int i = 0; i < cityCount; i++) initialPath[i] = i;

        double initialDist = tspEnergy.calculate(initialPath);
        System.out.printf("初始路径长度: %.2f%n", initialDist);

        // 模拟退火求解
        SimulatedAnnealing<int[]> sa = new SimulatedAnnealing<>(
            tspNeighbor, tspEnergy, exponentialCooling(0.995)
        );
        sa.setParameters(500.0, 10000, 5);

        SearchResult<int[]> result = sa.search(initialPath);
        System.out.printf("最优路径长度: %.2f%n", result.bestEnergy);
        System.out.printf("改进幅度: %.2f%%%n", (initialDist - result.bestEnergy) / initialDist * 100);

        // 展示最优路径前10个城市
        System.out.print("最优路径: ");
        for (int i = 0; i < Math.min(10, cityCount); i++) {
            System.out.print(result.bestSolution[i] + " -> ");
        }
        System.out.println("...");
    }

    public static void main(String[] args) {
        functionOptimizationDemo();
        tspDemo();
    }
}

复杂度分析

  • 单次迭代:生成邻居O(邻居构造复杂度),计算能量O(能量函数复杂度)
  • 总迭代次数:maxIterations × iterationsPerTemp
  • TSP案例:邻居生成O(n),能量计算O(n),总复杂度O(maxIterations × n)
  • 空间复杂度:O(n)存储当前解与最优解

参数调优指南

模拟退火的性能高度依赖参数选择:

  1. 初始温度:应使初期接受劣解的概率约为80%。可通过预实验统计初始阶段的平均ΔE来估算
  2. 降温系数:指数降温中α取0.95~0.995,迭代次数多则取大值
  3. 终止条件:温度低于阈值、连续多次无改进、或达到最大迭代次数
  4. 每个温度迭代次数:应与问题规模成正比,通常10~100次

与其他算法的对比

算法 特点 适用场景
模拟退火 单点搜索,概率跳出局部最优 连续/离散优化,解空间平滑
遗传算法 种群进化,多点并行搜索 组合优化,解空间复杂
禁忌搜索 记忆已访问解,避免循环 组合优化,邻域结构明确
粒子群 群体智能,收敛速度快 连续优化,实时性要求高

模拟退火以其简单优雅和广泛的适用性,成为解决NP难问题的首选启发式算法之一。掌握其核心原理与调参技巧,能为各类优化问题提供强有力的解决思路。