#!/usr/bin/env python3 # Simulated_Annealing.py """ Kapitel Metaheuristiken: Simulated Annealing an der Lackieranlage. Aufgabe: n Auftraege in eine Reihenfolge bringen. Zwischen zwei Auftraegen faellt eine Ruestzeit an - die Anlage muss gereinigt werden. Wie lange das dauert, haengt von BEIDEN Auftraegen ab: Ein Wechsel innerhalb derselben Produktfamilie kostet fast nichts, ein Wechsel von Dunkel nach Hell ist teuer, und zwischen manchen Familien ist der Reinigungsaufwand schlicht hoch. Gesucht ist die Reihenfolge mit der kleinsten Summe der Ruestzeiten. Das Programm zeigt drei Dinge: 1. Warum die Faustregel "immer der billigste naechste Auftrag" in eine Sackgasse laeuft - und wie weit man sie mit lokaler Suche verbessert. 2. Wie gross der Beitrag der Verschlechterungen WIRKLICH ist. Die Antwort faellt bescheidener aus als das Lehrbuch verspricht, und genau das ist der Grund, den Vergleich immer mitzurechnen. 3. Dass die Starttemperatur kein Beiwerk ist: zu heiss macht die Suche nutzlos, und man sieht es an einer einzigen Kennzahl kommen. Zum Budget: Der Abkuehlplan laeuft ueber eine feste ZUGZAHL, nicht ueber eine Sekundenzahl. Im Betrieb hat man zwar ein Zeitbudget - fuer einen VERGLEICH ist die Zugzahl aber das ehrlichere Mass: Sie ist auf jeder Maschine dieselbe, und die abgedruckten Zahlen unten lassen sich damit nachrechnen. Die gemessene Laufzeit steht trotzdem dabei. Benoetigt: numpy """ from __future__ import annotations import math import time import numpy as np ZUEGE = 400_000 # Zugbudget je Lauf (statt Sekunden: reproduzierbar) SAAT = 11 # feste Saat -> reproduzierbare Instanz # --- 1. Die Instanz --------------------------------------------------------- def erzeuge_ruestmatrix(n: int, saat: int = SAAT) -> np.ndarray: """Ruestzeiten in Minuten zwischen je zwei Auftraegen. Die Zahl der Produktfamilien waechst mit der Auftragszahl (eine Familie je zehn Auftraege). Sonst wuerde das Problem mit wachsendem n LEICHTER: Bei fester Familienzahl haette der Plan irgendwann nur noch ein paar grosse Bloecke, und jede Faustregel faende ihn. """ rng = np.random.default_rng(saat) familien = max(6, n // 10) # Reinigungsaufwand zwischen den Familien - unregelmaessig und asymmetrisch, # so wie im Betrieb: Von Klarlack auf Rot ist etwas anderes als umgekehrt. zwischen = rng.integers(8, 60, (familien, familien)) np.fill_diagonal(zwischen, 2) familie = rng.integers(0, familien, n) farbe = rng.integers(0, 10, n) # 0 = weiss ... 9 = schwarz matrix = np.zeros((n, n), dtype=np.int64) for i in range(n): for j in range(n): if i != j: # Dunkel -> hell kostet zusaetzlich: die Anlage muss heller # werden, das braucht mehr Spuelgaenge. 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: """Summe der Ruestzeiten einer Reihenfolge (offene Kette, keine Rundreise).""" return int(sum(matrix[reihe[k], reihe[k + 1]] for k in range(len(reihe) - 1))) def faustregel(matrix: np.ndarray) -> list[int]: """Immer der billigste noch offene Auftrag - 'naechster Nachbar'. So plant ein Meister von Hand, und es ist keine schlechte Regel. Ihr Fehler ist die Kurzsichtigkeit: Sie spart am Anfang und laesst die teuren Wechsel fuer das Ende uebrig, wo keine Wahl mehr bleibt. """ 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 # --- 2. Der Zug und seine Kostenaenderung ----------------------------------- def delta_verschieben(reihe: list[int], matrix: np.ndarray, von: int, nach: int): """Auftrag von Position 'von' an Position 'nach' versetzen. Liefert (Kostenaenderung, Restliste, Einfuegeposition) oder None, wenn der Zug nichts aendert. DAS IST DIE WICHTIGSTE FUNKTION DES PROGRAMMS. Sie berechnet die Kostenaenderung aus HOECHSTENS SECHS Matrixeintraegen, statt die Reihenfolge neu durchzusummieren. Der Unterschied ist nicht kosmetisch: Die volle Summe kostet O(n) je Zug, diese Rechnung O(1). Bei 500 Auftraegen sind das rund 500-mal mehr geprueste Zuege im selben Zeitbudget - und eine Metaheuristik lebt von der Zahl der Zuege. Herausgerissen wird der Auftrag aus zwei Kanten, seine Nachbarn ruecken zusammen (eine neue Kante). Eingefuegt wird er zwischen zwei andere Nachbarn, deren bisherige Kante dadurch verschwindet. """ 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: # die Nachbarn ruecken zusammen 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): # die aufgetrennte Kante faellt weg hinein -= matrix[rest[stelle - 1], rest[stelle]] return int(hinein - heraus), rest, stelle # --- 3. Die Suche ----------------------------------------------------------- class Verlauf: """Was ein Lauf ausser dem Ergebnis noch verraet.""" def __init__(self) -> None: self.schritte = 0 self.verschlechterungen_geprueft = 0 self.verschlechterungen_genommen = 0 self.schlechtester_zustand = 0 # wie weit die Suche abgedriftet ist self.sekunden = 0.0 # nur zur Information, nicht zum Vergleich @property def annahmequote(self) -> float: if not self.verschlechterungen_geprueft: return 0.0 return self.verschlechterungen_genommen / self.verschlechterungen_geprueft def suche(matrix: np.ndarray, zuege: int, start_temperatur: float, end_temperatur: float = 0.05, saat: int = 1, bergsteigen: bool = False) -> tuple[list[int], int, Verlauf]: """Simulated Annealing (oder reines Bergsteigen, wenn bergsteigen=True). Der einzige Unterschied zwischen beiden steht in der Annahmezeile: Das Bergsteigen nimmt nur Verbesserungen, SA nimmt Verschlechterungen mit einer Wahrscheinlichkeit an, die mit der Zeit gegen null geht. """ rng = np.random.default_rng(saat) n = len(matrix) reihe = faustregel(matrix) kosten = gesamtruestzeit(reihe, matrix) beste, beste_kosten = reihe[:], kosten verlauf = Verlauf() verlauf.schlechtester_zustand = kosten start = time.perf_counter() for zug in range(zuege): verlauf.schritte += 1 # Geometrischer Abkuehlplan ueber das Zugbudget. temperatur = start_temperatur * (end_temperatur / 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: annehmen = True else: verlauf.verschlechterungen_geprueft += 1 # Metropolis-Kriterium: je groesser die Verschlechterung und je # kaelter es ist, desto unwahrscheinlicher. annehmen = (not bergsteigen and rng.random() < math.exp(-aenderung / temperatur)) verlauf.verschlechterungen_genommen += annehmen if annehmen: reihe = rest[:stelle] + [reihe[von]] + rest[stelle:] kosten += aenderung verlauf.schlechtester_zustand = max(verlauf.schlechtester_zustand, kosten) if kosten < beste_kosten: beste, beste_kosten = reihe[:], kosten verlauf.sekunden = time.perf_counter() - start return beste, beste_kosten, verlauf if __name__ == "__main__": N = 200 matrix = erzeuge_ruestmatrix(N) start_reihe = faustregel(matrix) start_kosten = gesamtruestzeit(start_reihe, matrix) print("=" * 78) print(" SIMULATED ANNEALING AN DER LACKIERANLAGE") print("=" * 78) print(f"{N} Auftraege, {max(6, N // 10)} Produktfamilien, " f"Zugbudget {ZUEGE:,} je Lauf.\n") print(f"Faustregel 'billigster naechster Auftrag': " f"{start_kosten:,} Minuten Ruestzeit") # --- Wie gross sind die Zuege ueberhaupt? --------------------------- # Ohne diese Zahl kann man die Temperatur nicht waehlen: Das # Metropolis-Kriterium vergleicht die Verschlechterung MIT der Temperatur. rng = np.random.default_rng(0) stichprobe = [] for _ in range(5000): e = delta_verschieben(start_reihe, matrix, int(rng.integers(0, N)), int(rng.integers(0, N + 1))) if e: stichprobe.append(abs(e[0])) print(f"Typische Zuggroesse |Delta|: Median {np.median(stichprobe):.0f}, " f"90 %-Quantil {np.quantile(stichprobe, 0.9):.0f} Minuten") # --- Bergsteigen gegen Annealing ------------------------------------- print("\n" + "-" * 78) print("Nur Verbesserungen annehmen (Bergsteigen) gegen Annealing:\n") _, berg, berg_verlauf = suche(matrix, ZUEGE, 1.0, bergsteigen=True) _, gluehen, gluehen_verlauf = suche(matrix, ZUEGE, 1.0) print(f" {'Faustregel (Start)':<34} {start_kosten:>7,} Minuten") print(f" {'Bergsteigen':<34} {berg:>7,} Minuten " f"({(start_kosten - berg) / start_kosten * 100:5.1f} % besser)") print(f" {'Simulated Annealing':<34} {gluehen:>7,} Minuten " f"({(start_kosten - gluehen) / start_kosten * 100:5.1f} % besser)") print(f"\n Beide pruefen genau {berg_verlauf.schritte:,} Zuege " f"({berg_verlauf.sekunden:.1f} bzw. {gluehen_verlauf.sekunden:.1f} Sekunden).") print(f" Annealing nimmt davon {gluehen_verlauf.annahmequote * 100:.2f} % der " f"VERSCHLECHTERUNGEN an, Bergsteigen keine.") print("\n Bemerkenswert ist, wie klein der Abstand ist: Das simple") print(f" Bergsteigen holt hier bereits {(start_kosten - berg) / (start_kosten - gluehen) * 100:.0f} % dessen, was Annealing") print(" schafft. Wer eine Metaheuristik einsetzt, sollte diesen Vergleich") print(" IMMER mitrechnen - sonst schreibt man einer aufwendigen Methode") print(" gut, was schon die einfache geliefert haette.") # --- Die Temperatur ist kein Beiwerk --------------------------------- print("\n" + "-" * 78) print("Und warum die Starttemperatur ausgemessen gehoert:\n") print(f" {'Start-T':>8} {'Ruestzeit':>11} {'gg. Faustregel':>15} " f"{'angenommene':>12} {'schlechtester':>14}") print(f" {'':>8} {'':>11} {'':>15} {'Verschlecht.':>12} {'Zwischenwert':>14}") print(" " + "-" * 64) ergebnisse, verlaeufe = {}, {} for temperatur in (0.5, 1.0, 2.0, 4.0, 8.0): _, wert, verlauf = suche(matrix, ZUEGE, temperatur) ergebnisse[temperatur] = wert verlaeufe[temperatur] = verlauf gewinn = (start_kosten - wert) / start_kosten * 100 print(f" {temperatur:>8.1f} {wert:>11,} {gewinn:>13.1f} % " f"{verlauf.annahmequote * 100:>11.2f} % " f"{verlauf.schlechtester_zustand:>14,}") beste_temperatur = min(ergebnisse, key=ergebnisse.get) heiss = verlaeufe[8.0] print(f"\n Bester Wert bei T0 = {beste_temperatur}.") print("\n Die letzten beiden Spalten erklaeren, warum es nach oben kippt.") print(f" Bei T0 = 8 gehen nur {heiss.annahmequote * 100:.2f} % der Verschlechterungen") print(" durch - wenig genug, dass man es fuer harmlos halten koennte. Es") print(f" genuegt aber: Die Suche driftet bis auf {heiss.schlechtester_zustand:,} Minuten ab, das") print(f" {heiss.schlechtester_zustand / start_kosten:.1f}-fache der Startloesung, und findet im Zugbudget nicht") print(" mehr zurueck. Vom Ergebnis bleibt fast nichts uebrig: Die Suche hat") print(f" in {ZUEGE:,} Zuegen nie etwas GESEHEN, das mehr als") print(f" {(start_kosten - ergebnisse[8.0]) / start_kosten * 100:.1f} % unter der Startloesung lag.") print("\n Bemerkenswert ist die andere Richtung: Die besten Werte entstehen") print(" bei Annahmequoten nahe null. Die Lehrbuchregel 'anfangs sollen") print(" 20 bis 50 Prozent der Verschlechterungen durchgehen' fuehrt hier in") print(" die Irre - sie stammt aus Problemen mit sehr kleinen Zuggroessen.") print(" Bei Ruestzeiten in Minuten ist ein einziger schlechter Zug teuer.") print("\n Was stattdessen traegt: Den SCHLECHTESTEN Zwischenwert protokollieren") print(" und T0 so waehlen, dass er die Startloesung nicht wesentlich") print(" ueberschreitet. Diese eine Zahl haette hier auf Anhieb gezeigt, dass") print(" T0 = 8 unbrauchbar ist - ohne dass man das Endergebnis abwarten muss.") print("=" * 78)