// Copyright (C) 2026 Kiyotsugu Arai
// SPDX-License-Identifier: LGPL-3.0-or-later

// SimulatedAnnealing.hpp — 焼きなまし法 (Simulated Annealing)
#ifndef SANGI_OPTIMIZATION_SIMULATED_ANNEALING_HPP
#define SANGI_OPTIMIZATION_SIMULATED_ANNEALING_HPP

#include <vector>
#include <functional>
#include <random>
#include <cmath>

namespace sangi {

    // 焼きなまし法の結果
    template<typename T>
    struct SAResult {
        std::vector<T> bestSolution;
        T bestEnergy;
        size_t iterations;
    };

    // 焼きなまし法
    //   energy:      目的関数 (最小化)
    //   neighbor:    近傍生成関数
    //   x0:          初期解
    //   tempInit:    初期温度
    //   tempMin:     最低温度 (収束判定)
    //   coolingRate: 冷却率 (0 < coolingRate < 1)
    //   maxIter:     最大反復回数
    //   seed:        乱数シード
    template<typename T>
    SAResult<T> simulatedAnnealing(
        std::function<T(const std::vector<T>&)> energy,
        std::function<std::vector<T>(const std::vector<T>&, std::mt19937&)> neighbor,
        const std::vector<T>& x0,
        T tempInit = T(100), T tempMin = T(1e-8), T coolingRate = T(0.995),
        size_t maxIter = 100000, unsigned seed = 42)
    {
        std::mt19937 rng(seed);
        std::uniform_real_distribution<T> uniform(T(0), T(1));

        std::vector<T> current = x0;
        T currentEnergy = energy(current);

        std::vector<T> best = current;
        T bestEnergy = currentEnergy;

        T temp = tempInit;
        size_t iter = 0;

        for (; iter < maxIter && temp > tempMin; ++iter) {
            std::vector<T> candidate = neighbor(current, rng);
            T candidateEnergy = energy(candidate);
            T delta = candidateEnergy - currentEnergy;

            // Metropolis 判定: エネルギーが下がるか、確率 exp(-ΔE/T) で受理
            if (delta < T(0) || uniform(rng) < std::exp(-delta / temp)) {
                current = std::move(candidate);
                currentEnergy = candidateEnergy;

                if (currentEnergy < bestEnergy) {
                    best = current;
                    bestEnergy = currentEnergy;
                }
            }

            temp *= coolingRate;
        }

        return SAResult<T>{ std::move(best), bestEnergy, iter };
    }

} // namespace sangi

#endif // SANGI_OPTIMIZATION_SIMULATED_ANNEALING_HPP