# Pに属する問題の実装例
 
# 1. 2-SAT（多項式時間で解ける充足可能性問題）
from collections import defaultdict
 
class TwoSAT:
    """
    2-SAT問題の解法
    時間計算量: O(n + m)（含意グラフ上のSCC分解）
 
    2-SATはPに属する（多項式時間で解ける）
    一方、3-SATはNP完全
    """
 
    def __init__(self, n):
        """n: 変数の数"""
        self.n = n
        self.graph = defaultdict(list)
        self.reverse_graph = defaultdict(list)
 
    def _neg(self, x):
        """変数xの否定を返す"""
        return x + self.n if x < self.n else x - self.n
 
    def add_clause(self, x, y):
        """
        節 (x ∨ y) を追加
        含意: ¬x → y, ¬y → x
        """
        neg_x = self._neg(x)
        neg_y = self._neg(y)
        self.graph[neg_x].append(y)
        self.graph[neg_y].append(x)
        self.reverse_graph[y].append(neg_x)
        self.reverse_graph[x].append(neg_y)
 
    def solve(self):
        """
        SCC分解による2-SAT解法
 
        返り値: 充足可能な場合は各変数の真偽値のリスト、
                充足不可能な場合はNone
        """
        # Kosarajuのアルゴリズムでトポロジカル順序を求める
        visited = set()
        order = []
 
        def dfs1(v):
            visited.add(v)
            for u in self.graph[v]:
                if u not in visited:
                    dfs1(u)
            order.append(v)
 
        for v in range(2 * self.n):
            if v not in visited:
                dfs1(v)
 
        # 逆グラフでSCC分解
        comp = [-1] * (2 * self.n)
        comp_id = 0
 
        def dfs2(v, c):
            comp[v] = c
            for u in self.reverse_graph[v]:
                if comp[u] == -1:
                    dfs2(u, c)
 
        for v in reversed(order):
            if comp[v] == -1:
                dfs2(v, comp_id)
                comp_id += 1
 
        # 充足可能性チェック
        for i in range(self.n):
            if comp[i] == comp[i + self.n]:
                return None  # x と ¬x が同じSCC → 充足不可能
 
        # 解の構築
        result = []
        for i in range(self.n):
            result.append(comp[i] > comp[i + self.n])
 
        return result
 
 
# 使用例
sat = TwoSAT(3)  # 3変数: x₀, x₁, x₂
# (x₀ ∨ x₁) ∧ (¬x₀ ∨ x₂) ∧ (¬x₁ ∨ ¬x₂)
sat.add_clause(0, 1)      # x₀ ∨ x₁
sat.add_clause(3, 2)      # ¬x₀ ∨ x₂  (3 = ¬0)
sat.add_clause(4, 5)      # ¬x₁ ∨ ¬x₂  (4 = ¬1, 5 = ¬2)
 
solution = sat.solve()
print(f"解: {solution}")  # [True, True, False] など
 
 
# 2. 最大二部マッチング（Hopcroft-Karp）
class BipartiteMatching:
    """
    二部グラフの最大マッチング
    Hopcroft-Karpアルゴリズム: O(E√V)
    → Pに属する
    """
 
    def __init__(self, n, m):
        """n: 左頂点数、m: 右頂点数"""
        self.n = n
        self.m = m
        self.graph = defaultdict(list)
        self.match_left = [-1] * n
        self.match_right = [-1] * m
 
    def add_edge(self, u, v):
        """左のuから右のvへの辺を追加"""
        self.graph[u].append(v)
 
    def bfs(self):
        """BFSで増加パスの長さを計算"""
        from collections import deque
        dist = [-1] * self.n
        queue = deque()
 
        for u in range(self.n):
            if self.match_left[u] == -1:
                dist[u] = 0
                queue.append(u)
 
        found = False
        while queue:
            u = queue.popleft()
            for v in self.graph[u]:
                w = self.match_right[v]
                if w == -1:
                    found = True
                elif dist[w] == -1:
                    dist[w] = dist[u] + 1
                    queue.append(w)
 
        self.dist = dist
        return found
 
    def dfs(self, u):
        """DFSで増加パスを見つけてマッチングを更新"""
        for v in self.graph[u]:
            w = self.match_right[v]
            if w == -1 or (self.dist[w] == self.dist[u] + 1 and self.dfs(w)):
                self.match_left[u] = v
                self.match_right[v] = u
                return True
        self.dist[u] = -1
        return False
 
    def solve(self):
        """最大マッチングを求める"""
        matching = 0
        while self.bfs():
            for u in range(self.n):
                if self.match_left[u] == -1:
                    if self.dfs(u):
                        matching += 1
        return matching

# NP問題の検証器の実装例
 
class NPVerifier:
    """NP問題の検証器コレクション"""
 
    @staticmethod
    def verify_sat(formula, assignment):
        """
        SAT（充足可能性問題）の検証器
        formula: CNF形式の論理式 [[1, -2, 3], [-1, 3], ...]
                 正の数 = 変数、負の数 = 否定
        assignment: {変数番号: True/False}
 
        検証の計算量: O(n × m)（n=変数数、m=節数）
        → 多項式時間で検証可能 → NPに属する
        """
        for clause in formula:
            satisfied = False
            for literal in clause:
                var = abs(literal)
                value = assignment.get(var, False)
                if literal > 0 and value:
                    satisfied = True
                    break
                if literal < 0 and not value:
                    satisfied = True
                    break
            if not satisfied:
                return False
        return True
 
    @staticmethod
    def verify_hamiltonian_cycle(graph, cycle):
        """
        ハミルトン閉路の検証器
        graph: 隣接リスト {頂点: [隣接頂点, ...]}
        cycle: 頂点のリスト
 
        検証: O(V)
        """
        n = len(graph)
 
        # 全頂点を1回ずつ訪問しているか
        if len(cycle) != n:
            return False
        if set(cycle) != set(graph.keys()):
            return False
 
        # 各辺が存在するか
        for i in range(n):
            u = cycle[i]
            v = cycle[(i + 1) % n]
            if v not in graph[u]:
                return False
 
        return True
 
    @staticmethod
    def verify_graph_coloring(graph, coloring, k):
        """
        k-彩色の検証器
        graph: 隣接リスト
        coloring: {頂点: 色}
        k: 使用可能な色の数
 
        検証: O(V + E)
        """
        # 色の数がk以下か
        if len(set(coloring.values())) > k:
            return False
 
        # 全頂点に色が割り当てられているか
        if set(coloring.keys()) != set(graph.keys()):
            return False
 
        # 隣接する頂点が同じ色でないか
        for u in graph:
            for v in graph[u]:
                if coloring[u] == coloring[v]:
                    return False
 
        return True
 
    @staticmethod
    def verify_subset_sum(numbers, target, subset):
        """
        部分集合和問題の検証器
        numbers: 数の集合
        target: 目標の和
        subset: 選んだ要素のインデックス集合
 
        検証: O(n)
        """
        # subsetがnumbersの部分集合か
        if not all(0 <= i < len(numbers) for i in subset):
            return False
 
        # 和がtargetに等しいか
        return sum(numbers[i] for i in subset) == target
 
    @staticmethod
    def verify_clique(graph, clique, k):
        """
        k-クリークの検証器
        graph: 隣接リスト
        clique: 頂点の集合
        k: クリークの大きさ
 
        検証: O(k²)
        """
        if len(clique) != k:
            return False
 
        # 全ペアが辺で結ばれているか
        clique_list = list(clique)
        for i in range(len(clique_list)):
            for j in range(i + 1, len(clique_list)):
                u, v = clique_list[i], clique_list[j]
                if v not in graph.get(u, []):
                    return False
        return True
 
 
# 検証の実演
verifier = NPVerifier()
 
# SAT の検証
formula = [[1, -2, 3], [-1, 3], [2, -3]]  # (x₁ ∨ ¬x₂ ∨ x₃) ∧ ...
assignment = {1: True, 2: False, 3: True}
print(f"SAT検証: {verifier.verify_sat(formula, assignment)}")  # True
 
# 部分集合和の検証
numbers = [3, 7, 1, 8, 4]
target = 12
subset = {0, 2, 3}  # 3 + 1 + 8 = 12
print(f"部分集合和検証: {verifier.verify_subset_sum(numbers, target, subset)}")

# NP完全問題間の帰着の実装例
 
# 帰着1: 3-SAT → 独立集合
def reduce_3sat_to_independent_set(clauses, num_vars):
    """
    3-SATを独立集合問題に帰着する
 
    入力: 3-CNF式（各節が最大3リテラル）
    出力: グラフGと整数k
 
    (formula ∈ 3-SAT) ⟺ (G has independent set of size k)
    """
    # 各節の各リテラルに1つの頂点を作成
    vertices = []
    for i, clause in enumerate(clauses):
        for j, literal in enumerate(clause):
            vertices.append((i, j, literal))
 
    k = len(clauses)  # 独立集合のサイズ = 節の数
 
    # 辺の構築
    edges = []
    for idx1, v1 in enumerate(vertices):
        for idx2, v2 in enumerate(vertices):
            if idx1 >= idx2:
                continue
            # 同じ節内のリテラル間に辺を引く（1節から1つだけ選ぶ）
            if v1[0] == v2[0]:
                edges.append((idx1, idx2))
            # 矛盾するリテラル間に辺を引く（x と ¬x）
            elif v1[2] == -v2[2]:
                edges.append((idx1, idx2))
 
    return vertices, edges, k
 
 
# 帰着2: 独立集合 → 頂点被覆
def reduce_independent_set_to_vertex_cover(graph, n, k):
    """
    独立集合を頂点被覆に帰着する
 
    Gにサイズkの独立集合がある ⟺ Gにサイズ(n-k)の頂点被覆がある
 
    証明:
    SがGの独立集合 → V\Sが頂点被覆
    （SにはS内の2頂点を結ぶ辺がない → 全辺はV\Sの頂点を含む）
    """
    # グラフはそのまま、目標サイズを変換
    return graph, n - k
 
 
# 帰着3: 頂点被覆 → 集合被覆
def reduce_vertex_cover_to_set_cover(graph, n, k):
    """
    頂点被覆を集合被覆に帰着する
 
    全体集合U = 辺の集合
    各頂点vに対して、vに接続する辺の集合Sv
    k個以下のSvで全辺を被覆 ⟺ サイズk以下の頂点被覆が存在
    """
    # 全体集合 = 辺の集合
    universe = set()
    for u in graph:
        for v in graph[u]:
            edge = (min(u, v), max(u, v))
            universe.add(edge)
 
    # 各頂点の集合
    sets = {}
    for u in graph:
        sets[u] = set()
        for v in graph[u]:
            edge = (min(u, v), max(u, v))
            sets[u].add(edge)
 
    return universe, sets, k
 
 
# 帰着4: 部分集合和 → 分割問題
def reduce_subset_sum_to_partition(numbers, target):
    """
    部分集合和を分割問題に帰着する
 
    分割問題: 集合Sを2つの部分集合に分割して、
    両方の和が等しくなるか？
 
    帰着: S = numbers ∪ {|sum(numbers) - 2*target|}
    """
    total = sum(numbers)
    diff = abs(total - 2 * target)
    new_numbers = numbers + [diff]
 
    # new_numbersを和が等しい2つに分割可能
    # ⟺ numbersの部分集合の和がtargetに等しい
    return new_numbers

# 空間効率的なアルゴリズムの例
 
# Savitchの定理に基づく到達可能性判定
def reachability_savitch(graph, start, end, n):
    """
    Savitchの定理に基づく到達可能性判定
    空間計算量: O(log²n)（再帰の深さ × 各レベルの使用量）
 
    アイデア: 「startからendに2^i歩以内で到達できるか？」を
    中間地点midを全て試して再帰的に判定
    """
    def can_reach(u, v, steps):
        """uからvにstepsステップ以内で到達可能か？"""
        if steps == 0:
            return u == v
        if steps == 1:
            return u == v or v in graph.get(u, [])
 
        half = steps // 2
        # 全ての中間地点を試す
        for mid in range(n):
            if can_reach(u, mid, half) and can_reach(mid, v, steps - half):
                return True
        return False
 
    return can_reach(start, end, n)
 
 
# PSPACE完全: QBF（量化ブール式）
def solve_qbf(formula):
    """
    量化ブール式の評価
    PSPACE完全問題
 
    例: ∀x ∃y (x ∨ y) ∧ (¬x ∨ ¬y)
 
    空間効率的に解ける（多項式空間）が
    時間は指数的になりうる
    """
    if not formula.quantifiers:
        # 量化子がなくなったら、ブール式を評価
        return evaluate_boolean(formula.expression, formula.assignment)
 
    quantifier, variable = formula.quantifiers[0]
    remaining = formula.without_first_quantifier()
 
    if quantifier == 'forall':
        # ∀x: x=True と x=False の両方で真
        remaining.assignment[variable] = True
        result_true = solve_qbf(remaining)
        remaining.assignment[variable] = False
        result_false = solve_qbf(remaining)
        return result_true and result_false
 
    elif quantifier == 'exists':
        # ∃x: x=True または x=False のどちらかで真
        remaining.assignment[variable] = True
        result_true = solve_qbf(remaining)
        if result_true:
            return True
        remaining.assignment[variable] = False
        return solve_qbf(remaining)

# 確率的アルゴリズムの例
 
import random
 
# RP の例: Miller-Rabin素数判定
def miller_rabin(n, k=20):
    """
    Miller-Rabin素数判定
    RPに属する:
    - n が合成数 → 確率 ≥ 1 - 4^{-k} で「合成数」と答える
    - n が素数 → 確率 1 で「素数」と答える
    """
    if n < 2:
        return False
    if n < 4:
        return True
    if n % 2 == 0:
        return False
 
    # n - 1 = 2^r × d (dは奇数)
    r, d = 0, n - 1
    while d % 2 == 0:
        r += 1
        d //= 2
 
    for _ in range(k):
        a = random.randrange(2, n - 1)
        x = pow(a, d, n)
 
        if x == 1 or x == n - 1:
            continue
 
        for _ in range(r - 1):
            x = pow(x, 2, n)
            if x == n - 1:
                break
        else:
            return False  # 合成数
 
    return True  # おそらく素数
 
 
# BPP の例: ランダム化最小カット（Kargerのアルゴリズム）
def karger_min_cut(graph, n):
    """
    Kargerのランダム化最小カットアルゴリズム
    BPPに属する
 
    1回の実行で最小カットを見つける確率: ≥ 2/n²
    n²ln(n) 回繰り返せば高確率で正解
    """
    import copy
 
    # 辺の収縮
    vertices = list(range(n))
    edges = copy.deepcopy(graph)  # (u, v) のリスト
 
    # Union-Find
    parent = list(range(n))
 
    def find(x):
        while parent[x] != x:
            parent[x] = parent[parent[x]]
            x = parent[x]
        return x
 
    def union(x, y):
        px, py = find(x), find(y)
        parent[px] = py
 
    remaining = n
    while remaining > 2:
        # ランダムに辺を選んで収縮
        while True:
            u, v = random.choice(edges)
            if find(u) != find(v):
                break
 
        union(u, v)
        remaining -= 1
 
    # 残った2つのスーパー頂点間の辺を数える
    cut_size = 0
    for u, v in edges:
        if find(u) != find(v):
            cut_size += 1
 
    return cut_size
 
 
# ZPP の例: ランダム化クイックソート
def randomized_quicksort(arr):
    """
    ランダム化クイックソート
    ZPPに属する: 常に正しい結果、期待実行時間 O(n log n)
    最悪の場合 O(n²) だが確率は非常に低い
    """
    if len(arr) <= 1:
        return arr
 
    pivot = random.choice(arr)
    left = [x for x in arr if x < pivot]
    middle = [x for x in arr if x == pivot]
    right = [x for x in arr if x > pivot]
 
    return randomized_quicksort(left) + middle + randomized_quicksort(right)

# 1. ビット全探索（n ≤ 20）
def knapsack_brute_force(weights, values, capacity):
    """
    ナップサック問題のビット全探索
    O(2^n): n ≤ 20 程度まで実用的
    """
    n = len(weights)
    best_value = 0
    best_items = 0
 
    for mask in range(1 << n):
        total_weight = 0
        total_value = 0
        for i in range(n):
            if mask & (1 << i):
                total_weight += weights[i]
                total_value += values[i]
 
        if total_weight <= capacity and total_value > best_value:
            best_value = total_value
            best_items = mask
 
    return best_value, best_items
 
 
# 2. 動的計画法（指数時間だが高速化）
def tsp_dp(dist, n):
    """
    巡回セールスマン問題のDP解法（Held-Karp）
    O(n² × 2^n): ビット全探索の O(n! × n) より大幅に高速
    n ≤ 25 程度まで実用的
    """
    INF = float('inf')
    # dp[S][i] = 集合Sの都市を訪問し、最後にiにいる場合の最小コスト
    dp = [[INF] * n for _ in range(1 << n)]
    dp[1][0] = 0  # 都市0から開始
 
    for S in range(1 << n):
        for u in range(n):
            if dp[S][u] == INF:
                continue
            if not (S & (1 << u)):
                continue
            for v in range(n):
                if S & (1 << v):
                    continue
                new_S = S | (1 << v)
                dp[new_S][v] = min(dp[new_S][v], dp[S][u] + dist[u][v])
 
    # 全都市を訪問して0に戻る
    full = (1 << n) - 1
    result = min(dp[full][i] + dist[i][0] for i in range(1, n))
 
    return result
 
 
# 3. 分岐限定法（Branch and Bound）
def branch_and_bound_tsp(dist, n):
    """
    TSPの分岐限定法
    最悪 O(n!), 平均的にはDPより速いことが多い
    枝刈りの効果で大幅に探索を削減
    """
    import heapq
 
    INF = float('inf')
    best = INF
 
    def lower_bound(visited, current, cost):
        """下界の計算（簡易版: 最小辺の和）"""
        lb = cost
        for i in range(n):
            if i not in visited:
                min_edge = min(dist[i][j] for j in range(n) if j != i)
                lb += min_edge
        return lb
 
    # 優先度付きキュー: (下界, コスト, 現在地, 訪問済み)
    pq = [(lower_bound({0}, 0, 0), 0, 0, frozenset({0}))]
 
    while pq:
        lb, cost, current, visited = heapq.heappop(pq)
 
        if lb >= best:
            continue  # 枝刈り
 
        if len(visited) == n:
            total = cost + dist[current][0]
            best = min(best, total)
            continue
 
        for next_city in range(n):
            if next_city in visited:
                continue
            new_cost = cost + dist[current][next_city]
            new_visited = visited | {next_city}
            lb = lower_bound(new_visited, next_city, new_cost)
            if lb < best:
                heapq.heappush(pq, (lb, new_cost, next_city, new_visited))
 
    return best

# 近似アルゴリズムの実装例
 
# 1. 頂点被覆の2-近似
def vertex_cover_2approx(graph):
    """
    頂点被覆の2-近似アルゴリズム
 
    近似比: 2（最適解の2倍以内を保証）
    時間計算量: O(V + E)
 
    アルゴリズム:
    1. 辺を1つ選ぶ
    2. その両端点を被覆に追加
    3. 追加した頂点に接続する辺を全て削除
    4. 辺がなくなるまで繰り返す
    """
    cover = set()
    edges = set()
 
    for u in graph:
        for v in graph[u]:
            if u < v:
                edges.add((u, v))
 
    remaining_edges = set(edges)
 
    while remaining_edges:
        # 辺を1つ選ぶ
        u, v = next(iter(remaining_edges))
        cover.add(u)
        cover.add(v)
 
        # 両端点に接続する辺を削除
        remaining_edges = {
            (a, b) for (a, b) in remaining_edges
            if a != u and a != v and b != u and b != v
        }
 
    return cover
 
# 証明: 2-近似
# - 選ばれた辺の集合をMとする（|M|本）
# - Mの辺は互いに端点を共有しない（マッチング）
# - cover = 2|M| 頂点
# - 最適解は各辺から少なくとも1つの端点を含む → OPT ≥ |M|
# - cover = 2|M| ≤ 2 × OPT ∎
 
 
# 2. 集合被覆の貪欲近似（O(log n)-近似）
def greedy_set_cover(universe, sets):
    """
    集合被覆の貪欲アルゴリズム
 
    近似比: H(n) = Σ_{i=1}^{n} 1/i ≈ ln(n)
    → これが最善（P ≠ NP なら (1-ε)ln(n) 未満は不可能）
 
    時間計算量: O(|universe| × |sets|)
    """
    uncovered = set(universe)
    selected = []
    remaining_sets = list(sets.items())
 
    while uncovered:
        # 最も多くの未被覆要素をカバーする集合を選ぶ
        best_set = max(
            remaining_sets,
            key=lambda s: len(s[1] & uncovered)
        )
 
        if len(best_set[1] & uncovered) == 0:
            break  # これ以上カバーできない
 
        selected.append(best_set[0])
        uncovered -= best_set[1]
        remaining_sets.remove(best_set)
 
    return selected
 
 
# 3. TSPのChristofides近似（3/2-近似）
def christofides_tsp(dist, n):
    """
    Christofides のアルゴリズム（概念的実装）
 
    三角不等式を満たすTSPに対して 3/2-近似を保証
    これは多項式時間で達成される最良の近似比
 
    手順:
    1. 最小全域木Tを求める
    2. Tの奇数次頂点の最小重み完全マッチングMを求める
    3. T ∪ M でオイラー閉路を求める
    4. ショートカットしてハミルトン閉路にする
    """
    # 1. 最小全域木
    mst = compute_mst(dist, n)
 
    # 2. 奇数次頂点の最小完全マッチング
    odd_vertices = [v for v in range(n) if degree(mst, v) % 2 == 1]
    matching = min_weight_perfect_matching(dist, odd_vertices)
 
    # 3. MST + マッチングでオイラーグラフ
    euler_graph = combine(mst, matching)
    euler_tour = find_euler_tour(euler_graph)
 
    # 4. ショートカット（重複頂点を飛ばす）
    visited = set()
    hamiltonian = []
    for v in euler_tour:
        if v not in visited:
            visited.add(v)
            hamiltonian.append(v)
 
    return hamiltonian
 
# 近似比の証明:
# MST ≤ OPT（最適TSPから1辺を除くとMSTの上界）
# マッチング ≤ OPT/2（奇数次頂点上の最適TSPの半分以下）
# よって: Christofides ≤ MST + Matching ≤ OPT + OPT/2 = 3/2 × OPT
 
 
# 4. ナップサック問題のFPTAS（完全多項式時間近似スキーム）
def knapsack_fptas(weights, values, capacity, epsilon):
    """
    ナップサック問題のFPTAS
 
    任意のε > 0 に対して、(1-ε)-近似を保証
    時間計算量: O(n² / ε)
 
    FPTASが存在 → 強いNP困難ではない
    """
    n = len(weights)
    v_max = max(values)
 
    # スケーリング
    K = epsilon * v_max / n
    scaled_values = [int(v / K) for v in values]
    V_total = sum(scaled_values)
 
    # スケールされた値でDP
    INF = float('inf')
    # dp[v] = 価値合計がちょうどvとなる最小重量
    dp = [INF] * (V_total + 1)
    dp[0] = 0
 
    for i in range(n):
        for v in range(V_total, scaled_values[i] - 1, -1):
            if dp[v - scaled_values[i]] + weights[i] < dp[v]:
                dp[v] = dp[v - scaled_values[i]] + weights[i]
 
    # 容量以内で最大の価値を見つける
    best_v = 0
    for v in range(V_total + 1):
        if dp[v] <= capacity:
            best_v = v
 
    return best_v * K  # 元のスケールに戻す

# ヒューリスティックの実装例
 
import random
import math
 
# 1. 焼きなまし法（Simulated Annealing）
def simulated_annealing_tsp(dist, n, initial_temp=10000,
                             cooling_rate=0.9995, min_temp=1e-8):
    """
    TSPの焼きなまし法
 
    保証: なし（ただし実用上は非常に良い解が得られる）
    実行時間: ユーザーが指定（反復回数 or 温度で制御）
    """
    # 初期解: ランダムな順列
    current = list(range(n))
    random.shuffle(current)
 
    def tour_cost(tour):
        return sum(dist[tour[i]][tour[(i+1) % n]] for i in range(n))
 
    current_cost = tour_cost(current)
    best = current[:]
    best_cost = current_cost
    temp = initial_temp
 
    while temp > min_temp:
        # 近傍: 2-opt（2辺を交換）
        i = random.randint(0, n - 2)
        j = random.randint(i + 1, n - 1)
        neighbor = current[:i] + current[i:j+1][::-1] + current[j+1:]
        neighbor_cost = tour_cost(neighbor)
 
        # 受理判定
        delta = neighbor_cost - current_cost
        if delta < 0 or random.random() < math.exp(-delta / temp):
            current = neighbor
            current_cost = neighbor_cost
 
            if current_cost < best_cost:
                best = current[:]
                best_cost = current_cost
 
        temp *= cooling_rate
 
    return best, best_cost
 
 
# 2. 遺伝的アルゴリズム（Genetic Algorithm）
def genetic_algorithm_tsp(dist, n, pop_size=100, generations=1000,
                          mutation_rate=0.02):
    """
    TSPの遺伝的アルゴリズム
 
    構成要素:
    - 個体: 都市の順列（染色体）
    - 適応度: 経路長の逆数
    - 交叉: 順序交叉（OX）
    - 突然変異: 2-opt
    - 選択: トーナメント選択
    """
    def tour_cost(tour):
        return sum(dist[tour[i]][tour[(i+1) % n]] for i in range(n))
 
    def tournament_select(population, fitnesses, tournament_size=5):
        indices = random.sample(range(len(population)), tournament_size)
        best_idx = max(indices, key=lambda i: fitnesses[i])
        return population[best_idx]
 
    def order_crossover(parent1, parent2):
        """順序交叉（OX）"""
        start = random.randint(0, n - 2)
        end = random.randint(start + 1, n - 1)
 
        child = [-1] * n
        child[start:end+1] = parent1[start:end+1]
 
        fill_pos = (end + 1) % n
        parent2_pos = (end + 1) % n
        while -1 in child:
            if parent2[parent2_pos] not in child:
                child[fill_pos] = parent2[parent2_pos]
                fill_pos = (fill_pos + 1) % n
            parent2_pos = (parent2_pos + 1) % n
 
        return child
 
    def mutate(tour):
        """2-opt突然変異"""
        if random.random() < mutation_rate:
            i = random.randint(0, n - 2)
            j = random.randint(i + 1, n - 1)
            tour[i:j+1] = reversed(tour[i:j+1])
        return tour
 
    # 初期集団
    population = [random.sample(range(n), n) for _ in range(pop_size)]
 
    best_ever = None
    best_cost_ever = float('inf')
 
    for gen in range(generations):
        # 適応度計算
        costs = [tour_cost(ind) for ind in population]
        fitnesses = [1.0 / c for c in costs]
 
        # 最良個体の更新
        best_idx = min(range(pop_size), key=lambda i: costs[i])
        if costs[best_idx] < best_cost_ever:
            best_cost_ever = costs[best_idx]
            best_ever = population[best_idx][:]
 
        # 次世代の生成
        new_population = [best_ever[:]]  # エリート保存
 
        while len(new_population) < pop_size:
            p1 = tournament_select(population, fitnesses)
            p2 = tournament_select(population, fitnesses)
            child = order_crossover(p1, p2)
            child = mutate(child)
            new_population.append(child)
 
        population = new_population
 
    return best_ever, best_cost_ever
 
 
# 3. タブー探索
def tabu_search_tsp(dist, n, max_iter=10000, tabu_size=20):
    """
    TSPのタブー探索
 
    直近のtabu_size回の移動を禁止（タブーリスト）
    → 局所最適解から脱出可能
    """
    current = list(range(n))
    random.shuffle(current)
 
    def tour_cost(tour):
        return sum(dist[tour[i]][tour[(i+1) % n]] for i in range(n))
 
    current_cost = tour_cost(current)
    best = current[:]
    best_cost = current_cost
    tabu_list = []
 
    for _ in range(max_iter):
        # 全ての2-opt近傍を列挙
        best_neighbor = None
        best_neighbor_cost = float('inf')
        best_move = None
 
        for i in range(n - 1):
            for j in range(i + 1, n):
                move = (i, j)
                neighbor = current[:i] + current[i:j+1][::-1] + current[j+1:]
                nc = tour_cost(neighbor)
 
                # タブーでない、または最良解を更新する場合
                if (move not in tabu_list or nc < best_cost):
                    if nc < best_neighbor_cost:
                        best_neighbor = neighbor
                        best_neighbor_cost = nc
                        best_move = move
 
        if best_neighbor is None:
            break
 
        current = best_neighbor
        current_cost = best_neighbor_cost
 
        # タブーリスト更新
        tabu_list.append(best_move)
        if len(tabu_list) > tabu_size:
            tabu_list.pop(0)
 
        if current_cost < best_cost:
            best = current[:]
            best_cost = current_cost
 
    return best, best_cost

# FPTアルゴリズムの例
 
# 頂点被覆のFPTアルゴリズム（分岐限定法）
def vertex_cover_fpt(graph, k):
    """
    頂点被覆のFPTアルゴリズム
    時間計算量: O(2^k × (V + E))
 
    アイデア:
    辺(u,v)を1つ選ぶ → uを被覆に入れるか、vを入れるかの2択
    → 深さkの二分木を探索
    """
    def solve(edges, cover, remaining_k):
        # 辺がなくなったら成功
        if not edges:
            return cover
 
        # 予算がなくなったら失敗
        if remaining_k == 0:
            return None
 
        # 辺を1つ選ぶ
        u, v = next(iter(edges))
 
        # 選択1: uを被覆に入れる
        new_edges_u = {(a, b) for (a, b) in edges if a != u and b != u}
        result = solve(new_edges_u, cover | {u}, remaining_k - 1)
        if result is not None:
            return result
 
        # 選択2: vを被覆に入れる
        new_edges_v = {(a, b) for (a, b) in edges if a != v and b != v}
        result = solve(new_edges_v, cover | {v}, remaining_k - 1)
        if result is not None:
            return result
 
        return None  # サイズk以下の頂点被覆は存在しない
 
    # 辺の集合を構築
    edges = set()
    for u in graph:
        for v in graph[u]:
            if u < v:
                edges.add((u, v))
 
    return solve(edges, set(), k)
 
 
# Color Coding: k-パス問題のFPTアルゴリズム
def color_coding_k_path(graph, n, k):
    """
    k-パス問題: グラフに長さkのパスが存在するか？
    FPTアルゴリズム: O(2^k × E) （ランダム化）
 
    アイデア:
    1. 各頂点にランダムに k 色を割り当て
    2. 全ての色が使われるカラフルなパスをDPで見つける
    3. 成功確率 ≥ e^{-k} → O(e^k) 回繰り返せば高確率で見つかる
    """
    for trial in range(int(2.72 ** k) * 2):  # e^k × 2 回試行
        # ランダムな彩色
        colors = [random.randint(0, k - 1) for _ in range(n)]
 
        # DP: dp[v][S] = 頂点vで終わり、色集合Sを使うカラフルなパスが存在するか
        dp = [[False] * (1 << k) for _ in range(n)]
 
        for v in range(n):
            dp[v][1 << colors[v]] = True
 
        for S in range(1, 1 << k):
            for v in range(n):
                if not dp[v][S]:
                    continue
                for u in graph.get(v, []):
                    c = colors[u]
                    if not (S & (1 << c)):  # 新しい色
                        dp[u][S | (1 << c)] = True
 
        # 全色を使うパスが見つかったか？
        full_set = (1 << k) - 1
        for v in range(n):
            if dp[v][full_set]:
                return True
 
    return False

"""
演習: 貪欲法による集合被覆の近似アルゴリズムを実装し、近似比を実測せよ。
"""
 
import random
 
def generate_set_cover_instance(n, m, density=0.3):
    """集合被覆のランダムインスタンス生成"""
    universe = set(range(n))
    sets = {}
    for i in range(m):
        size = max(1, int(n * density * random.random()))
        sets[i] = set(random.sample(list(universe), min(size, n)))
 
    # 全要素が被覆可能であることを保証
    for elem in universe:
        random_set = random.choice(list(sets.keys()))
        sets[random_set].add(elem)
 
    return universe, sets
 
 
def optimal_set_cover(universe, sets):
    """厳密解（小規模用、指数時間）"""
    n = len(sets)
    keys = list(sets.keys())
    best = None
 
    for mask in range(1, 1 << n):
        covered = set()
        selected = []
        for i in range(n):
            if mask & (1 << i):
                covered |= sets[keys[i]]
                selected.append(keys[i])
        if covered >= universe:
            if best is None or len(selected) < len(best):
                best = selected
 
    return best
 
 
# 実験
for n in [10, 15, 20]:
    total_ratio = 0
    trials = 50
    for _ in range(trials):
        universe, sets = generate_set_cover_instance(n, n * 2)
        greedy = greedy_set_cover(universe, sets)
        optimal = optimal_set_cover(universe, sets)
        ratio = len(greedy) / len(optimal)
        total_ratio += ratio
 
    avg_ratio = total_ratio / trials
    print(f"n={n}: 平均近似比 = {avg_ratio:.3f} (理論上界: {math.log(n):.3f})")

"""
演習: DPLLアルゴリズム（SATソルバーの基礎）を実装せよ。
"""
 
def dpll(clauses, assignment=None):
    """
    DPLLアルゴリズム — SAT問題の標準的な解法
 
    最悪の場合 O(2^n) だが、実用上は非常に効率的
    現代のSATソルバー（MiniSat, Z3等）の基礎
    """
    if assignment is None:
        assignment = {}
 
    # 充足チェック
    clauses = simplify(clauses, assignment)
 
    # 空の節がある → 充足不可能
    if any(len(c) == 0 for c in clauses):
        return None
 
    # 全ての節が除去された → 充足可能
    if len(clauses) == 0:
        return assignment
 
    # ユニット伝播（Unit Propagation）
    unit_clauses = [c for c in clauses if len(c) == 1]
    while unit_clauses:
        literal = unit_clauses[0][0]
        var = abs(literal)
        value = literal > 0
        assignment[var] = value
        clauses = simplify(clauses, {var: value})
 
        if any(len(c) == 0 for c in clauses):
            return None
        if len(clauses) == 0:
            return assignment
 
        unit_clauses = [c for c in clauses if len(c) == 1]
 
    # 純リテラル除去（Pure Literal Elimination）
    all_literals = set()
    for clause in clauses:
        for lit in clause:
            all_literals.add(lit)
 
    for lit in list(all_literals):
        if -lit not in all_literals:
            var = abs(lit)
            assignment[var] = (lit > 0)
            clauses = simplify(clauses, {var: lit > 0})
 
    if len(clauses) == 0:
        return assignment
 
    # 分岐（Branching）
    var = abs(clauses[0][0])
 
    # True を試す
    result = dpll(clauses, {**assignment, var: True})
    if result is not None:
        return result
 
    # False を試す
    result = dpll(clauses, {**assignment, var: False})
    return result
 
 
def simplify(clauses, assignment):
    """割り当てに基づいて節を簡約化"""
    new_clauses = []
    for clause in clauses:
        new_clause = []
        satisfied = False
        for lit in clause:
            var = abs(lit)
            if var in assignment:
                if (lit > 0) == assignment[var]:
                    satisfied = True
                    break
                # リテラルがFalseの場合は節から除去
            else:
                new_clause.append(lit)
        if not satisfied:
            new_clauses.append(new_clause)
    return new_clauses