#!/usr/bin/env python3 # Metaheuristik_vs_Exakt.py """ Kapitel Metaheuristiken: Ab wann lohnt sich eine Metaheuristik? Die Antwort ist eine Zahl. Die verbreitete Regel lautet: "Kleine Probleme exakt, grosse mit Metaheuristiken." Sie ist richtig - aber voellig nutzlos, solange niemand sagt, wo 'gross' anfaengt. Dieses Programm misst es fuer die Ruestzeitenaufgabe aus Simulated_Annealing.py nach, und das Ergebnis ueberrascht in beide Richtungen: * Bei 60 Auftraegen loest CP-SAT BEWEISBAR optimal, in wenigen Sekunden. Die Metaheuristik ist hier schlicht die schlechtere Wahl - sie liefert ein schwaecheres Ergebnis und weiss nicht einmal, wie schwach. * Bei 200 Auftraegen liegt CP-SAT ohne Optimalitaetsbeweis vorne. * Erst bei 500 Auftraegen kippt es: Dort gibt CP-SAT im Zeitlimit eine Loesung aus, die SCHLECHTER ist als die Faustregel eines Meisters. Der Punkt des Kapitels steht in der letzten Spalte: Die Schranke. Nur der exakte Solver sagt, wie gut die Loesung sein KANN. Eine Metaheuristik allein liefert eine Zahl ohne Massstab - man weiss nie, ob man 2 % oder 40 % neben dem Optimum liegt. ZUR REPRODUZIERBARKEIT: Die Spalte 'Annealing' ist auf jeder Maschine gleich - sie laeuft ueber ein festes ZUGbudget. Die CP-SAT-Spalten sind es nicht: Ein Zeitlimit ist ein Wanduhr-Limit, und wie weit der Solver in 30 Sekunden kommt, haengt von der Maschine ab. CP-SAT kennt zwar ein deterministisches Zeitmass (max_deterministic_time), das aber in einer Groessenordnung rechnet, die fuer ein Buchbeispiel unbrauchbar langsam ist. Die Zahlen unten sind also hardwareabhaengig - die AUSSAGE der Tabelle ist es nicht, sie ist auf verschiedenen Groessen und in mehreren Laeufen stabil. WICHTIG: Hier wird nur ortools importiert, nicht highspy (Kapitel Oekosystem). Benoetigt: numpy, ortools """ from __future__ import annotations import time import numpy as np from ortools.sat.python import cp_model ZUEGE = 400_000 # Zugbudget der Metaheuristik (wie Simulated_Annealing.py) ZEITLIMIT = 30.0 # Zeitbudget des exakten Solvers GROESSEN = (60, 200, 500) # --- Instanz, Faustregel und Suche: identisch zu Simulated_Annealing.py ----- def erzeuge_ruestmatrix(n: int, saat: int = 11) -> np.ndarray: rng = np.random.default_rng(saat) familien = max(6, n // 10) zwischen = rng.integers(8, 60, (familien, familien)) np.fill_diagonal(zwischen, 2) familie = rng.integers(0, familien, n) farbe = rng.integers(0, 10, n) matrix = np.zeros((n, n), dtype=np.int64) for i in range(n): for j in range(n): if i != j: matrix[i, j] = (zwischen[familie[i], familie[j]] + 3 * max(0, farbe[i] - farbe[j])) return matrix def gesamtruestzeit(reihe: list[int], matrix: np.ndarray) -> int: return int(sum(matrix[reihe[k], reihe[k + 1]] for k in range(len(reihe) - 1))) def faustregel(matrix: np.ndarray) -> list[int]: n = len(matrix) offen = set(range(1, n)) reihe = [0] while offen: naechster = min(offen, key=lambda j: matrix[reihe[-1], j]) reihe.append(naechster) offen.discard(naechster) return reihe def delta_verschieben(reihe, matrix, von, nach): n = len(reihe) if von == nach or nach == von + 1: return None heraus = 0 if von > 0: heraus += matrix[reihe[von - 1], reihe[von]] if von < n - 1: heraus += matrix[reihe[von], reihe[von + 1]] if 0 < von < n - 1: heraus -= matrix[reihe[von - 1], reihe[von + 1]] rest = reihe[:von] + reihe[von + 1:] stelle = nach if nach < von else nach - 1 hinein = 0 if stelle > 0: hinein += matrix[rest[stelle - 1], reihe[von]] if stelle < len(rest): hinein += matrix[reihe[von], rest[stelle]] if 0 < stelle < len(rest): hinein -= matrix[rest[stelle - 1], rest[stelle]] return int(hinein - heraus), rest, stelle def annealing(matrix: np.ndarray, zuege: int = ZUEGE, start_temperatur: float = 1.0, saat: int = 1) -> tuple[int, float]: """Dieselbe Suche wie in Simulated_Annealing.py, hier nur als Vergleichswert.""" import math rng = np.random.default_rng(saat) n = len(matrix) reihe = faustregel(matrix) kosten = gesamtruestzeit(reihe, matrix) beste_kosten = kosten t0 = time.perf_counter() for zug in range(zuege): temperatur = start_temperatur * (0.05 / start_temperatur) ** (zug / zuege) von = int(rng.integers(0, n)) nach = int(rng.integers(0, n + 1)) ergebnis = delta_verschieben(reihe, matrix, von, nach) if ergebnis is None: continue aenderung, rest, stelle = ergebnis if aenderung <= 0 or rng.random() < math.exp(-aenderung / temperatur): reihe = rest[:stelle] + [reihe[von]] + rest[stelle:] kosten += aenderung beste_kosten = min(beste_kosten, kosten) return beste_kosten, time.perf_counter() - t0 # --- Der exakte Weg: CP-SAT mit AddCircuit --------------------------------- def exakt(matrix: np.ndarray, zeitlimit: float = ZEITLIMIT): """Reihenfolge als Rundreise mit einem kostenlosen Hilfsknoten. AddCircuit verlangt einen geschlossenen Kreis. Unsere Aufgabe ist aber eine offene Kette: Der erste Auftrag hat keinen Vorgaenger, der letzte keinen Nachfolger. Der uebliche Kniff ist ein zusaetzlicher Hilfsknoten n, dessen Kanten nichts kosten. Der Kreis laeuft dann Hilfsknoten -> Kette -> Hilfsknoten, und die Kosten des Kreises sind genau die der Kette. AddCircuit ist der Grund, warum CP-SAT hier so lange mithaelt: Die Nebenbedingung sorgt selbst dafuer, dass keine Teilkreise entstehen. Wer dasselbe als MILP formuliert, braucht dafuer entweder exponentiell viele Ungleichungen oder die schwache MTZ-Formulierung (Kapitel Graphen). """ n = len(matrix) modell = cp_model.CpModel() kanten, ziel = [], [] for i in range(n + 1): for j in range(n + 1): if i == j: continue aktiv = modell.NewBoolVar(f"kante_{i}_{j}") kanten.append((i, j, aktiv)) if i < n and j < n: # Hilfsknoten kostet nichts ziel.append(int(matrix[i, j]) * aktiv) modell.AddCircuit(kanten) modell.Minimize(sum(ziel)) loeser = cp_model.CpSolver() loeser.parameters.max_time_in_seconds = zeitlimit loeser.parameters.num_workers = 8 loeser.parameters.random_seed = 1 t0 = time.perf_counter() status = loeser.Solve(modell) dauer = time.perf_counter() - t0 return (loeser.StatusName(status), int(loeser.ObjectiveValue()), int(loeser.BestObjectiveBound()), dauer) if __name__ == "__main__": print("=" * 88) print(" AB WANN LOHNT SICH DIE METAHEURISTIK?") print("=" * 88) print(f"Dieselbe Aufgabe in drei Groessen. Metaheuristik: {ZUEGE:,} Zuege. " f"Exakt: {ZEITLIMIT:.0f} s Zeitlimit.\n") print(f"{'Auftraege':>10} {'Faustregel':>11} {'Annealing':>10} " f"{'CP-SAT':>9} {'Status':>10} {'Schranke':>9} {'CP-Zeit':>9}") print("-" * 88) zeilen = [] for n in GROESSEN: matrix = erzeuge_ruestmatrix(n) start = gesamtruestzeit(faustregel(matrix), matrix) heuristisch, _ = annealing(matrix) status, wert, schranke, dauer = exakt(matrix) zeilen.append((n, start, heuristisch, wert, status, schranke, dauer)) print(f"{n:>10} {start:>11,} {heuristisch:>10,} {wert:>9,} " f"{status:>10} {schranke:>9,} {dauer:>8.1f}s") print("-" * 88) print("\nWas in jeder Zeile steht:\n") for n, start, heur, wert, status, schranke, _ in zeilen: if status == "OPTIMAL": urteil = (f"CP-SAT beweist das Optimum. Annealing liegt " f"{(heur - wert) / wert * 100:.1f} % daneben - und haette") zweite = " das ohne den exakten Lauf nicht erfahren." elif wert <= heur: urteil = (f"CP-SAT liegt {(heur - wert) / heur * 100:.1f} % vor Annealing, " f"ohne Beweis (Gap {(wert - schranke) / wert * 100:.1f} %).") zweite = " Auch ohne Optimalitaetsbeweis ist es die bessere Wahl." else: urteil = (f"UMSCHLAGPUNKT: CP-SAT ist {(wert - heur) / heur * 100:.1f} % " f"SCHLECHTER als Annealing") zweite = (f" und sogar {(wert - start) / start * 100:.1f} % schlechter " f"als die Faustregel des Meisters.") print(f" {n:>4} Auftraege: {urteil}") print(zweite) n, start, heur, wert, status, schranke, _ = zeilen[-1] print("\n" + "=" * 88) print(" DIE SPALTE, DIE MAN NICHT WEGLASSEN DARF") print("=" * 88) print("Die Schranke ist der eigentliche Ertrag des exakten Solvers - auch dann,") print("wenn seine LOESUNG unbrauchbar ist.") print() print(f"Bei {n} Auftraegen liefert Annealing {heur:,} Minuten. Klingt gut: " f"{(start - heur) / start * 100:.1f} %") print("besser als die Faustregel. Die Schranke aus demselben CP-SAT-Lauf, dessen") print(f"Loesung wir gerade verworfen haben, sagt aber: Unter {schranke:,} Minuten geht") print(f"es nicht. Zwischen {schranke:,} und {heur:,} liegen {(heur - schranke) / heur * 100:.0f} % - so viel Luft") print("kann noch nach oben sein.") print() print("Das ist der Unterschied zwischen 'besser als vorher' und 'gut'. Eine") print("Metaheuristik allein kann nur das Erste. Deshalb gehoert auch dann ein") print("exakter Lauf dazu, wenn man am Ende die heuristische Loesung einsetzt:") print("nicht wegen seiner Loesung, sondern wegen seiner Schranke.") print("=" * 88)