m1une's library

This documentation is automatically generated by online-judge-tools/verification-helper

View on GitHub

:heavy_check_mark: verify/heuristic/simulated_annealing.test.cpp

Depends on

Code

#define PROBLEM "https://judge.yosupo.jp/problem/aplusb"

#include <algorithm>
#include <cassert>
#include <cmath>
#include <iostream>
#include <limits>

#include "../../heuristic/all.hpp"
#include "../../utilities/random.hpp"

using m1une::heuristic::AnnealingCooling;
using m1une::heuristic::AnnealingObjective;
using m1une::heuristic::SimulatedAnnealing;

bool close(double first, double second) {
    return std::abs(first - second) <=
           1e-12 * std::max(1.0, std::abs(second));
}

void test_temperature() {
    SimulatedAnnealing exponential(100.0, 1.0);
    assert(close(exponential.temperature(0.0), 100.0));
    assert(close(exponential.temperature(0.5), 10.0));
    assert(close(exponential.temperature(1.0), 1.0));

    SimulatedAnnealing linear(100.0, 0.0,
                              AnnealingObjective::maximize,
                              AnnealingCooling::linear);
    assert(close(linear.temperature(0.0), 100.0));
    assert(close(linear.temperature(0.25), 75.0));
    assert(close(linear.temperature(1.0), 0.0));
}

void test_maximization() {
    SimulatedAnnealing annealing(10.0, 10.0);
    assert(annealing.acceptance_probability(3, 4, 0.5) == 1.0);
    assert(annealing.acceptance_probability(3, 3, 0.5) == 1.0);
    assert(close(annealing.acceptance_probability(3, -7, 0.5),
                 std::exp(-1.0)));
    assert(annealing.accept(3, 4, 0.5, 0.999999));
    assert(annealing.accept(3, -7, 0.5, 0.3));
    assert(!annealing.accept(3, -7, 0.5, 0.4));
}

void test_minimization() {
    SimulatedAnnealing annealing(2.0, 2.0,
                                 AnnealingObjective::minimize);
    assert(annealing.acceptance_probability(5, 4, 0.0) == 1.0);
    assert(close(annealing.acceptance_probability(5, 7, 0.0),
                 std::exp(-1.0)));
    assert(annealing.accept_delta(-1.0, 0.0, 0.999999));
    assert(!annealing.accept_delta(2.0, 0.0, 0.4));
}

void test_zero_temperature_and_large_scores() {
    SimulatedAnnealing greedy(1.0, 0.0,
                              AnnealingObjective::maximize,
                              AnnealingCooling::linear);
    assert(greedy.acceptance_probability_delta(-1.0, 1.0) == 0.0);
    assert(!greedy.accept_delta(-1.0, 1.0, 0.0));
    assert(greedy.accept_delta(0.0, 1.0, 0.999999));

    long long low = std::numeric_limits<long long>::min();
    long long high = std::numeric_limits<long long>::max();
    assert(greedy.acceptance_probability(low, high, 1.0) == 1.0);
    assert(greedy.acceptance_probability(high, low, 1.0) == 0.0);
}

void test_randomized_against_formula() {
    m1une::utilities::Random random(0x51a7edULL);
    for (AnnealingObjective objective : {AnnealingObjective::minimize,
                                         AnnealingObjective::maximize}) {
        SimulatedAnnealing annealing(30.0, 0.03, objective);
        for (int trial = 0; trial < 10000; trial++) {
            long long current = random.uniform(-1000000000, 1000000000);
            long long candidate = random.uniform(-1000000000, 1000000000);
            double progress = random.real();
            double random01 = random.real();

            long double delta = static_cast<long double>(candidate) -
                                static_cast<long double>(current);
            long double improvement =
                objective == AnnealingObjective::maximize ? delta : -delta;
            double expected = 1.0;
            if (improvement < 0.0L) {
                expected = std::exp(static_cast<double>(
                    improvement / annealing.temperature(progress)));
            }
            assert(close(annealing.acceptance_probability(
                             current, candidate, progress),
                         expected));
            assert(annealing.accept(current, candidate, progress, random01) ==
                   (random01 < expected));
        }
    }
}

int main() {
    test_temperature();
    test_maximization();
    test_minimization();
    test_zero_temperature_and_large_scores();
    test_randomized_against_formula();

    long long a, b;
    std::cin >> a >> b;
    std::cout << a + b << '\n';
}
#line 1 "verify/heuristic/simulated_annealing.test.cpp"
#define PROBLEM "https://judge.yosupo.jp/problem/aplusb"

#include <algorithm>
#include <cassert>
#include <cmath>
#include <iostream>
#include <limits>

#line 1 "heuristic/all.hpp"



#line 1 "heuristic/beam_search.hpp"



#line 6 "heuristic/beam_search.hpp"
#include <concepts>
#include <cstddef>
#include <functional>
#include <type_traits>
#include <utility>
#include <vector>

#line 1 "heuristic/objective.hpp"



namespace m1une {
namespace heuristic {

enum class Objective {
    minimize,
    maximize,
};

template <class Score>
bool better_score(const Score& first, const Score& second,
                  Objective objective) {
    if (objective == Objective::maximize) return second < first;
    return first < second;
}

}  // namespace heuristic
}  // namespace m1une


#line 14 "heuristic/beam_search.hpp"

namespace m1une {
namespace heuristic {

template <class State, class Score>
struct BeamSearchResult {
    State state;
    Score score;
    int depth;
    std::size_t expanded_states;
    std::size_t generated_states;
};

namespace beam_search_detail {

template <class State, class Score>
struct Node {
    State state;
    Score score;
    std::size_t order;
};

template <class State, class Score>
struct BetterNode {
    Objective objective;

    bool operator()(const Node<State, Score>& first,
                    const Node<State, Score>& second) const {
        if (better_score(first.score, second.score, objective)) return true;
        if (better_score(second.score, first.score, objective)) return false;
        return first.order < second.order;
    }
};

}  // namespace beam_search_detail

// expand(state, next_depth) may return a range of children. For allocation-free
// generation, expand(state, next_depth, emit) may instead call emit(child).
// evaluate(state) returns its score. The best beam_width states are retained at
// every depth, and the best state in the last non-empty layer is returned.
template <class State, class Expand, class Evaluate>
auto beam_search(State initial_state, int depth_limit, int beam_width,
                 Expand expand, Evaluate evaluate,
                 Objective objective = Objective::maximize) {
    assert(0 <= depth_limit);
    assert(0 < beam_width);

    using Score = std::remove_cvref_t<
        std::invoke_result_t<Evaluate&, const State&>>;
    using Node = beam_search_detail::Node<State, Score>;
    using Better = beam_search_detail::BetterNode<State, Score>;

    Score initial_score = std::invoke(evaluate, initial_state);
    std::vector<Node> beam;
    beam.push_back(Node{std::move(initial_state),
                        std::move(initial_score), 0});

    std::size_t expanded_states = 0;
    std::size_t generated_states = 0;
    int reached_depth = 0;
    if (depth_limit < 0 || beam_width <= 0) depth_limit = 0;

    Better better{objective};
    for (int next_depth = 1; next_depth <= depth_limit; next_depth++) {
        std::vector<Node> candidates;
        candidates.reserve(static_cast<std::size_t>(beam_width));
        std::size_t order = 0;

        for (const Node& node : beam) {
            expanded_states++;
            auto emit = [&](auto&& candidate_state) {
                using Candidate = decltype(candidate_state);
                static_assert(std::is_constructible_v<State, Candidate>);
                State state(std::forward<Candidate>(candidate_state));
                Score candidate_score = std::invoke(evaluate, state);
                Node candidate{std::move(state), std::move(candidate_score),
                               order++};
                generated_states++;
                if (int(candidates.size()) < beam_width) {
                    candidates.push_back(std::move(candidate));
                    std::push_heap(candidates.begin(), candidates.end(), better);
                } else if (better(candidate, candidates.front())) {
                    std::pop_heap(candidates.begin(), candidates.end(), better);
                    candidates.back() = std::move(candidate);
                    std::push_heap(candidates.begin(), candidates.end(), better);
                }
            };
            if constexpr (std::invocable<Expand&, const State&, int>) {
                auto next_states =
                    std::invoke(expand, node.state, next_depth);
                for (auto& candidate_state : next_states) {
                    emit(std::move(candidate_state));
                }
            } else if constexpr (std::invocable<Expand&, const State&>) {
                auto next_states = std::invoke(expand, node.state);
                for (auto& candidate_state : next_states) {
                    emit(std::move(candidate_state));
                }
            } else {
                std::invoke(expand, node.state, next_depth, emit);
            }
        }

        if (candidates.empty()) break;
        beam = std::move(candidates);
        reached_depth = next_depth;
    }

    int best = 0;
    for (int index = 1; index < int(beam.size()); index++) {
        if (better(beam[index], beam[best])) best = index;
    }
    return BeamSearchResult<State, Score>{
        std::move(beam[best].state), std::move(beam[best].score),
        reached_depth, expanded_states, generated_states};
}

}  // namespace heuristic
}  // namespace m1une


#line 1 "heuristic/hill_climbing.hpp"



#line 5 "heuristic/hill_climbing.hpp"

#line 7 "heuristic/hill_climbing.hpp"

namespace m1une {
namespace heuristic {

using HillClimbingObjective = Objective;

class HillClimbing {
   private:
    Objective _objective;
    bool _accept_equal;

   public:
    explicit HillClimbing(Objective objective = Objective::maximize,
                          bool accept_equal = false)
        : _objective(objective), _accept_equal(accept_equal) {}

    bool accept_delta(long double candidate_minus_current) const {
        if (_objective == Objective::maximize) {
            return _accept_equal ? 0.0L <= candidate_minus_current
                                 : 0.0L < candidate_minus_current;
        }
        return _accept_equal ? candidate_minus_current <= 0.0L
                             : candidate_minus_current < 0.0L;
    }

    template <std::convertible_to<long double> CurrentScore,
              std::convertible_to<long double> CandidateScore>
    bool accept(CurrentScore current_score,
                CandidateScore candidate_score) const {
        long double delta = static_cast<long double>(candidate_score) -
                            static_cast<long double>(current_score);
        return accept_delta(delta);
    }
};

}  // namespace heuristic
}  // namespace m1une


#line 1 "heuristic/simulated_annealing.hpp"



#line 8 "heuristic/simulated_annealing.hpp"

#line 10 "heuristic/simulated_annealing.hpp"

namespace m1une {
namespace heuristic {

using AnnealingObjective = Objective;

enum class AnnealingCooling {
    linear,
    exponential,
};

class SimulatedAnnealing {
   private:
    double _start_temperature;
    double _end_temperature;
    AnnealingObjective _objective;
    AnnealingCooling _cooling;

    long double directed_delta(long double candidate_minus_current) const {
        if (_objective == AnnealingObjective::maximize) {
            return candidate_minus_current;
        }
        return -candidate_minus_current;
    }

   public:
    SimulatedAnnealing(
        double start_temperature, double end_temperature,
        AnnealingObjective objective = AnnealingObjective::maximize,
        AnnealingCooling cooling = AnnealingCooling::exponential)
        : _start_temperature(start_temperature),
          _end_temperature(end_temperature),
          _objective(objective),
          _cooling(cooling) {
        assert(std::isfinite(start_temperature));
        assert(std::isfinite(end_temperature));
        assert(0.0 <= end_temperature);
        assert(end_temperature <= start_temperature);
        assert(cooling != AnnealingCooling::exponential ||
               0.0 < end_temperature);
    }

    double temperature(double progress) const {
        assert(std::isfinite(progress));
        assert(0.0 <= progress && progress <= 1.0);
        progress = std::clamp(progress, 0.0, 1.0);
        if (_cooling == AnnealingCooling::linear) {
            return _start_temperature +
                   (_end_temperature - _start_temperature) * progress;
        }
        return _start_temperature *
               std::pow(_end_temperature / _start_temperature, progress);
    }

    double acceptance_probability_delta(
        long double candidate_minus_current, double progress) const {
        long double improvement = directed_delta(candidate_minus_current);
        if (0.0L <= improvement) return 1.0;
        double current_temperature = temperature(progress);
        if (current_temperature == 0.0) return 0.0;
        return std::exp(static_cast<double>(
            improvement / static_cast<long double>(current_temperature)));
    }

    bool accept_delta(long double candidate_minus_current, double progress,
                      double random01) const {
        assert(std::isfinite(random01));
        assert(0.0 <= random01 && random01 < 1.0);
        return random01 <
               acceptance_probability_delta(candidate_minus_current, progress);
    }

    template <std::convertible_to<long double> CurrentScore,
              std::convertible_to<long double> CandidateScore>
    double acceptance_probability(CurrentScore current_score,
                                  CandidateScore candidate_score,
                                  double progress) const {
        long double delta = static_cast<long double>(candidate_score) -
                            static_cast<long double>(current_score);
        return acceptance_probability_delta(delta, progress);
    }

    template <std::convertible_to<long double> CurrentScore,
              std::convertible_to<long double> CandidateScore>
    bool accept(CurrentScore current_score, CandidateScore candidate_score,
                double progress, double random01) const {
        long double delta = static_cast<long double>(candidate_score) -
                            static_cast<long double>(current_score);
        return accept_delta(delta, progress, random01);
    }
};

}  // namespace heuristic
}  // namespace m1une


#line 8 "heuristic/all.hpp"


#line 1 "utilities/random.hpp"



#line 6 "utilities/random.hpp"
#include <chrono>
#line 8 "utilities/random.hpp"
#include <cstdint>
#line 10 "utilities/random.hpp"
#include <numeric>
#include <queue>
#include <random>
#include <string>
#include <string_view>
#include <tuple>
#line 17 "utilities/random.hpp"
#include <unordered_set>
#line 20 "utilities/random.hpp"

namespace m1une {
namespace utilities {

struct RandomGraphOptions {
    bool directed = false;
    bool allow_self_loops = false;
    bool allow_parallel_edges = false;
};

struct Random {
   private:
    std::mt19937_64 _engine;

    static unsigned long long chrono_seed() {
        return static_cast<unsigned long long>(
            std::chrono::steady_clock::now().time_since_epoch().count());
    }

    static std::uint64_t graph_edge_count(int vertex_count,
                                          const RandomGraphOptions& options) {
        std::uint64_t n = static_cast<unsigned int>(vertex_count);
        if (options.directed) {
            return options.allow_self_loops ? n * n : n * (n - 1);
        }
        return options.allow_self_loops ? n * (n + 1) / 2 : n * (n - 1) / 2;
    }

    static std::pair<int, int> decode_graph_edge(
        std::uint64_t index, int vertex_count,
        const RandomGraphOptions& options) {
        std::uint64_t n = static_cast<unsigned int>(vertex_count);
        if (options.directed) {
            std::uint64_t width = options.allow_self_loops ? n : n - 1;
            int from = int(index / width);
            int offset = int(index % width);
            int to = options.allow_self_loops || offset < from ? offset : offset + 1;
            return {from, to};
        }

        auto prefix = [&](std::uint64_t vertex) {
            if (options.allow_self_loops) {
                return vertex * (2 * n - vertex + 1) / 2;
            }
            return vertex * (2 * n - vertex - 1) / 2;
        };
        std::uint64_t low = 0;
        std::uint64_t high = n;
        while (low + 1 < high) {
            std::uint64_t middle = (low + high) / 2;
            if (prefix(middle) <= index) {
                low = middle;
            } else {
                high = middle;
            }
        }
        int from = int(low);
        int to = from + int(index - prefix(low)) +
                 (options.allow_self_loops ? 0 : 1);
        return {from, to};
    }

   public:
    Random() : _engine(chrono_seed()) {}
    explicit Random(unsigned long long seed) : _engine(seed) {}

    void seed(unsigned long long value) {
        _engine.seed(value);
    }

    std::mt19937_64& engine() {
        return _engine;
    }

    unsigned long long operator()() {
        return _engine();
    }

    long long uniform(long long l, long long r) {
        return std::uniform_int_distribution<long long>(l, r)(_engine);
    }

    unsigned long long uniform_unsigned(unsigned long long l, unsigned long long r) {
        return std::uniform_int_distribution<unsigned long long>(l, r)(_engine);
    }

    double real(double l = 0.0, double r = 1.0) {
        return std::uniform_real_distribution<double>(l, r)(_engine);
    }

    template <std::integral T>
    requires(!std::same_as<std::remove_cv_t<T>, bool>)
    std::vector<T> sequence(int size, T lower, T upper) {
        assert(0 <= size);
        assert(lower <= upper);
        if (size < 0 || upper < lower) return {};
        std::vector<T> result(size);
        if constexpr (std::signed_integral<T>) {
            std::uniform_int_distribution<long long> distribution(
                static_cast<long long>(lower), static_cast<long long>(upper));
            for (T& value : result) value = static_cast<T>(distribution(_engine));
        } else {
            std::uniform_int_distribution<unsigned long long> distribution(
                static_cast<unsigned long long>(lower),
                static_cast<unsigned long long>(upper));
            for (T& value : result) value = static_cast<T>(distribution(_engine));
        }
        return result;
    }

    std::string string(
        int length,
        std::string_view alphabet = "abcdefghijklmnopqrstuvwxyz") {
        assert(0 <= length);
        assert(length == 0 || !alphabet.empty());
        if (length < 0 || (0 < length && alphabet.empty())) return {};
        std::string result(length, '\0');
        for (char& character : result) {
            character = alphabet[uniform(0, int(alphabet.size()) - 1)];
        }
        return result;
    }

    std::vector<int> permutation(int size, int first = 0) {
        assert(0 <= size);
        if (size < 0) return {};
        std::vector<int> result(size);
        std::iota(result.begin(), result.end(), first);
        shuffle(result);
        return result;
    }

    // Returns the edges of a uniformly random labeled tree on [0, size).
    std::vector<std::pair<int, int>> tree(int size) {
        assert(0 <= size);
        if (size <= 1) return {};

        std::vector<int> prufer = sequence(size - 2, 0, size - 1);
        std::vector<int> degree(size, 1);
        for (int vertex : prufer) degree[vertex]++;
        std::priority_queue<int, std::vector<int>, std::greater<int>> leaves;
        for (int vertex = 0; vertex < size; vertex++) {
            if (degree[vertex] == 1) leaves.push(vertex);
        }

        std::vector<std::pair<int, int>> edges;
        edges.reserve(size - 1);
        for (int vertex : prufer) {
            int leaf = leaves.top();
            leaves.pop();
            edges.emplace_back(leaf, vertex);
            if (--degree[vertex] == 1) leaves.push(vertex);
        }
        int first = leaves.top();
        leaves.pop();
        edges.emplace_back(first, leaves.top());

        shuffle(edges);
        for (auto& [from, to] : edges) {
            if (uniform(0, 1)) std::swap(from, to);
        }
        return edges;
    }

    // Returns m random edges on [0, vertex_count). By default the result is
    // a simple undirected graph without self-loops.
    std::vector<std::pair<int, int>> graph(
        int vertex_count, int edge_count,
        RandomGraphOptions options = {}) {
        assert(0 <= vertex_count);
        assert(0 <= edge_count);
        if (vertex_count < 0 || edge_count < 0) return {};
        if (edge_count == 0) return {};
        assert(0 < vertex_count);
        if (vertex_count == 0) return {};
        if (!options.allow_self_loops) {
            assert(2 <= vertex_count || edge_count == 0);
            if (vertex_count < 2) return {};
        }

        std::vector<std::pair<int, int>> edges;
        edges.reserve(edge_count);
        if (options.allow_parallel_edges) {
            for (int edge = 0; edge < edge_count; edge++) {
                int from = int(uniform(0, vertex_count - 1));
                int to;
                if (options.allow_self_loops) {
                    to = int(uniform(0, vertex_count - 1));
                } else {
                    to = int(uniform(0, vertex_count - 2));
                    if (from <= to) to++;
                }
                if (!options.directed && to < from) std::swap(from, to);
                edges.emplace_back(from, to);
            }
            return edges;
        }

        std::uint64_t maximum = graph_edge_count(vertex_count, options);
        assert(static_cast<std::uint64_t>(edge_count) <= maximum);
        if (maximum < static_cast<std::uint64_t>(edge_count)) return {};

        std::unordered_set<std::uint64_t> selected;
        selected.reserve(static_cast<std::size_t>(edge_count) * 2 + 1);
        std::vector<std::uint64_t> indices;
        indices.reserve(edge_count);
        for (std::uint64_t current = maximum - edge_count;
             current < maximum; current++) {
            std::uint64_t candidate = uniform_unsigned(0, current);
            if (selected.contains(candidate)) candidate = current;
            selected.insert(candidate);
            indices.push_back(candidate);
        }
        for (std::uint64_t index : indices) {
            edges.push_back(decode_graph_edge(index, vertex_count, options));
        }
        return edges;
    }

    std::vector<std::pair<int, int>> directed_graph(
        int vertex_count, int edge_count,
        bool allow_self_loops = false) {
        RandomGraphOptions options;
        options.allow_self_loops = allow_self_loops;
        return directed_graph(vertex_count, edge_count, options);
    }

    std::vector<std::pair<int, int>> directed_graph(
        int vertex_count, int edge_count, RandomGraphOptions options) {
        options.directed = true;
        return graph(vertex_count, edge_count, options);
    }

    // Returns a directed acyclic graph. Vertices are randomly permuted before
    // every sampled edge is directed forward in that topological order.
    std::vector<std::pair<int, int>> dag(
        int vertex_count, int edge_count,
        RandomGraphOptions options = {}) {
        options.directed = false;
        options.allow_self_loops = false;
        std::vector<std::pair<int, int>> edges =
            graph(vertex_count, edge_count, options);
        std::vector<int> order = permutation(vertex_count);
        for (auto& [from, to] : edges) {
            from = order[from];
            to = order[to];
        }
        return edges;
    }

    template <std::integral Weight>
    requires(!std::same_as<std::remove_cv_t<Weight>, bool>)
    std::vector<std::tuple<int, int, Weight>> weighted_tree(
        int size, Weight lower, Weight upper) {
        std::vector<std::pair<int, int>> edges = tree(size);
        std::vector<Weight> weights = sequence(int(edges.size()), lower, upper);
        std::vector<std::tuple<int, int, Weight>> result;
        result.reserve(edges.size());
        for (int index = 0; index < int(edges.size()); index++) {
            result.emplace_back(edges[index].first, edges[index].second,
                                weights[index]);
        }
        return result;
    }

    template <std::integral Weight>
    requires(!std::same_as<std::remove_cv_t<Weight>, bool>)
    std::vector<std::tuple<int, int, Weight>> weighted_graph(
        int vertex_count, int edge_count, Weight lower, Weight upper,
        RandomGraphOptions options = {}) {
        std::vector<std::pair<int, int>> edges =
            graph(vertex_count, edge_count, options);
        std::vector<Weight> weights = sequence(int(edges.size()), lower, upper);
        std::vector<std::tuple<int, int, Weight>> result;
        result.reserve(edges.size());
        for (int index = 0; index < int(edges.size()); index++) {
            result.emplace_back(edges[index].first, edges[index].second,
                                weights[index]);
        }
        return result;
    }

    template <std::integral Weight>
    requires(!std::same_as<std::remove_cv_t<Weight>, bool>)
    std::vector<std::tuple<int, int, Weight>> weighted_directed_graph(
        int vertex_count, int edge_count, Weight lower, Weight upper,
        bool allow_self_loops = false) {
        RandomGraphOptions options;
        options.allow_self_loops = allow_self_loops;
        return weighted_directed_graph(vertex_count, edge_count, lower, upper,
                                       options);
    }

    template <std::integral Weight>
    requires(!std::same_as<std::remove_cv_t<Weight>, bool>)
    std::vector<std::tuple<int, int, Weight>> weighted_directed_graph(
        int vertex_count, int edge_count, Weight lower, Weight upper,
        RandomGraphOptions options) {
        options.directed = true;
        return weighted_graph(vertex_count, edge_count, lower, upper, options);
    }

    template <std::integral Weight>
    requires(!std::same_as<std::remove_cv_t<Weight>, bool>)
    std::vector<std::tuple<int, int, Weight>> weighted_dag(
        int vertex_count, int edge_count, Weight lower, Weight upper,
        RandomGraphOptions options = {}) {
        std::vector<std::pair<int, int>> edges =
            dag(vertex_count, edge_count, options);
        std::vector<Weight> weights = sequence(int(edges.size()), lower, upper);
        std::vector<std::tuple<int, int, Weight>> result;
        result.reserve(edges.size());
        for (int index = 0; index < int(edges.size()); index++) {
            result.emplace_back(edges[index].first, edges[index].second,
                                weights[index]);
        }
        return result;
    }

    template <typename T>
    void shuffle(std::vector<T>& v) {
        std::shuffle(v.begin(), v.end(), _engine);
    }

    template <typename Iterator>
    void shuffle(Iterator first, Iterator last) {
        std::shuffle(first, last, _engine);
    }

    template <typename T>
    const T& choice(const std::vector<T>& v) {
        return v[uniform(0, static_cast<long long>(v.size()) - 1)];
    }
};

}  // namespace utilities
}  // namespace m1une


#line 11 "verify/heuristic/simulated_annealing.test.cpp"

using m1une::heuristic::AnnealingCooling;
using m1une::heuristic::AnnealingObjective;
using m1une::heuristic::SimulatedAnnealing;

bool close(double first, double second) {
    return std::abs(first - second) <=
           1e-12 * std::max(1.0, std::abs(second));
}

void test_temperature() {
    SimulatedAnnealing exponential(100.0, 1.0);
    assert(close(exponential.temperature(0.0), 100.0));
    assert(close(exponential.temperature(0.5), 10.0));
    assert(close(exponential.temperature(1.0), 1.0));

    SimulatedAnnealing linear(100.0, 0.0,
                              AnnealingObjective::maximize,
                              AnnealingCooling::linear);
    assert(close(linear.temperature(0.0), 100.0));
    assert(close(linear.temperature(0.25), 75.0));
    assert(close(linear.temperature(1.0), 0.0));
}

void test_maximization() {
    SimulatedAnnealing annealing(10.0, 10.0);
    assert(annealing.acceptance_probability(3, 4, 0.5) == 1.0);
    assert(annealing.acceptance_probability(3, 3, 0.5) == 1.0);
    assert(close(annealing.acceptance_probability(3, -7, 0.5),
                 std::exp(-1.0)));
    assert(annealing.accept(3, 4, 0.5, 0.999999));
    assert(annealing.accept(3, -7, 0.5, 0.3));
    assert(!annealing.accept(3, -7, 0.5, 0.4));
}

void test_minimization() {
    SimulatedAnnealing annealing(2.0, 2.0,
                                 AnnealingObjective::minimize);
    assert(annealing.acceptance_probability(5, 4, 0.0) == 1.0);
    assert(close(annealing.acceptance_probability(5, 7, 0.0),
                 std::exp(-1.0)));
    assert(annealing.accept_delta(-1.0, 0.0, 0.999999));
    assert(!annealing.accept_delta(2.0, 0.0, 0.4));
}

void test_zero_temperature_and_large_scores() {
    SimulatedAnnealing greedy(1.0, 0.0,
                              AnnealingObjective::maximize,
                              AnnealingCooling::linear);
    assert(greedy.acceptance_probability_delta(-1.0, 1.0) == 0.0);
    assert(!greedy.accept_delta(-1.0, 1.0, 0.0));
    assert(greedy.accept_delta(0.0, 1.0, 0.999999));

    long long low = std::numeric_limits<long long>::min();
    long long high = std::numeric_limits<long long>::max();
    assert(greedy.acceptance_probability(low, high, 1.0) == 1.0);
    assert(greedy.acceptance_probability(high, low, 1.0) == 0.0);
}

void test_randomized_against_formula() {
    m1une::utilities::Random random(0x51a7edULL);
    for (AnnealingObjective objective : {AnnealingObjective::minimize,
                                         AnnealingObjective::maximize}) {
        SimulatedAnnealing annealing(30.0, 0.03, objective);
        for (int trial = 0; trial < 10000; trial++) {
            long long current = random.uniform(-1000000000, 1000000000);
            long long candidate = random.uniform(-1000000000, 1000000000);
            double progress = random.real();
            double random01 = random.real();

            long double delta = static_cast<long double>(candidate) -
                                static_cast<long double>(current);
            long double improvement =
                objective == AnnealingObjective::maximize ? delta : -delta;
            double expected = 1.0;
            if (improvement < 0.0L) {
                expected = std::exp(static_cast<double>(
                    improvement / annealing.temperature(progress)));
            }
            assert(close(annealing.acceptance_probability(
                             current, candidate, progress),
                         expected));
            assert(annealing.accept(current, candidate, progress, random01) ==
                   (random01 < expected));
        }
    }
}

int main() {
    test_temperature();
    test_maximization();
    test_minimization();
    test_zero_temperature_and_large_scores();
    test_randomized_against_formula();

    long long a, b;
    std::cin >> a >> b;
    std::cout << a + b << '\n';
}
Back to top page