#!/usr/bin/env python3 # erzeuge_polyeder.py """ Erzeugt die beiden Polyeder-Diagramme zum Fundament- und zum LP-Kapitel: bilder_04/kap_fundament_polyeder.svg der zulaessige Bereich bilder_04/kap_lp_simplex_ecken.svg der Eckenlauf des Simplex Beide zeigen DASSELBE Modell - das Bot-Allokationsproblem, das sich durch die ersten Kapitel zieht: max 150 x1 + 250 x2 u.d.N. 2 x1 + 5 x2 <= 40 (vCPU) 4 x1 + 6 x2 <= 60 (RAM) x1 <= 8 (Marktlimit) x1, x2 >= 0 Warum das wichtig ist: Die Vorgaengerbilder zeigten ein voellig anderes LP (2x1+3x2 <= 24, 3x1+2x2 <= 24, Optimum bei (4,8; 4,8)) - unmittelbar neben einer Handrechnung, die mit dem Bot-Modell rechnet und (7,5; 5) herausbekommt. Ein Leser musste annehmen, beides gehoere zusammen. Die Ecken werden NICHT hingeschrieben, sondern nach demselben Verfahren bestimmt wie in Visualisierung_Loesungsraum.py: Schnittpunkte je zweier Begrenzungsgeraden, die zulaessig sind. Der Simplex-Weg im zweiten Bild ist der Weg der Handrechnung im LP-Kapitel: (0,0) -> (0,8) -> (7,5; 5), also drei der fuenf Ecken. Aufruf (aus dem Repository-Wurzelverzeichnis): python3 bilder_04/erzeuge_polyeder.py Benoetigt: numpy, matplotlib """ from __future__ import annotations import itertools import os import matplotlib matplotlib.use("Agg") import matplotlib.pyplot as plt import numpy as np # Reproduzierbare SVG-Ausgabe (siehe erzeuge_titelseite.py) plt.rcParams["svg.hashsalt"] = "or-mit-python-v04" BASIS = os.path.dirname(os.path.dirname(os.path.abspath(__file__))) BILDER = os.path.join(BASIS, "bilder_04") INDIGO, CYAN, AMBER, GRUEN = "#4338ca", "#0891b2", "#b45309", "#059669" VIOLETT = "#7c3aed" # Niveaulinien und Zielfunktionsgradient # Modell wortgleich zu Visualisierung_Loesungsraum.py C = np.array([150.0, 250.0]) A = np.array([[2.0, 5.0], [4.0, 6.0], [1.0, 0.0]]) B = np.array([40.0, 60.0, 8.0]) NAMEN = ["vCPU: $2x_1+5x_2 \\leq 40$", "RAM: $4x_1+6x_2 \\leq 60$", "Marktlimit: $x_1 \\leq 8$"] FARBEN = [INDIGO, CYAN, AMBER] # Der Weg der Handrechnung im LP-Kapitel: Start im Ursprung, dann x2 aufnehmen # (Pivotspalte -250), dann x1 (Pivotspalte -50). SIMPLEXWEG = [(0.0, 0.0), (0.0, 8.0), (7.5, 5.0)] def ist_zulaessig(punkt: np.ndarray, tol: float = 1e-7) -> bool: return bool(np.all(A @ punkt <= B + tol) and np.all(punkt >= -tol)) def berechne_ecken() -> np.ndarray: """Ecken = zulaessige Schnittpunkte je zweier Begrenzungsgeraden. Wortgleich zu Visualisierung_Loesungsraum.py - die Koordinaten im Bild stammen also aus derselben Rechnung wie die im Buch abgedruckte Tabelle, nicht aus einer zweiten Quelle. """ geraden = [(A[i], B[i]) for i in range(len(B))] geraden.append((np.array([1.0, 0.0]), 0.0)) # x1 = 0 geraden.append((np.array([0.0, 1.0]), 0.0)) # x2 = 0 ecken: list[np.ndarray] = [] for (n1, d1), (n2, d2) in itertools.combinations(geraden, 2): matrix = np.array([n1, n2]) if abs(np.linalg.det(matrix)) < 1e-9: continue punkt = np.linalg.solve(matrix, np.array([d1, d2])) if ist_zulaessig(punkt) and not any(np.allclose(punkt, e) for e in ecken): ecken.append(punkt) return np.array(ecken) # Beschriftungsposition je Ecke, in Punkten relativ zum Eckpunkt. Von Hand # gesetzt und nicht gerechnet: (7,5; 5) und (8; 4,67) liegen eine halbe Einheit # auseinander, und jede automatische Regel legt ihre Beschriftungen uebereinander. # Das ist Layout, keine Aussage - die Koordinaten selbst bleiben gerechnet. VERSATZ = { (0.0, 0.0): (14, 10), (0.0, 8.0): (14, -2), (8.0, 0.0): (-10, 14), (7.5, 5.0): (-12, 16), (8.0, 4.67): (16, -4), } # Im Simplex-Bild laeuft der Pfad durch zwei dieser Positionen - dort andere. VERSATZ_SIMPLEX = { (0.0, 8.0): (16, 12), (7.5, 5.0): (14, 16), } def versatz(punkt: np.ndarray, mit_pfad: bool = False) -> tuple: schluessel = (round(float(punkt[0]), 2), round(float(punkt[1]), 2)) if mit_pfad and schluessel in VERSATZ_SIMPLEX: return VERSATZ_SIMPLEX[schluessel] return VERSATZ.get(schluessel, (12, 8)) def punkt_text(punkt: np.ndarray) -> str: """Koordinaten deutsch und auf zwei Stellen - '(8; 4,67)' statt '(8; 4.66667)'.""" teile = [] for wert in punkt: text = f"{wert:.2f}".rstrip("0").rstrip(".") teile.append(text.replace(".", ",")) return f"({teile[0]}; {teile[1]})" def euro(wert: float) -> str: return f"{wert:,.0f} €".replace(",", ".") def sortiere_umlaufend(ecken: np.ndarray) -> np.ndarray: """Ecken gegen den Uhrzeigersinn ordnen, damit das Polygon nicht verknotet.""" mitte = ecken.mean(axis=0) winkel = np.arctan2(ecken[:, 1] - mitte[1], ecken[:, 0] - mitte[0]) return ecken[np.argsort(winkel)] def zeichne_rahmen(achse, ecken: np.ndarray) -> None: """Restriktionsgeraden, zulaessiger Bereich und Achsenbeschriftung.""" rand = sortiere_umlaufend(ecken) achse.fill(rand[:, 0], rand[:, 1], color=INDIGO, alpha=0.10, zorder=1) achse.plot(np.append(rand[:, 0], rand[0, 0]), np.append(rand[:, 1], rand[0, 1]), color=INDIGO, linewidth=1.6, zorder=2) x = np.linspace(0, 11, 200) for zeile, grenze, name, farbe in zip(A, B, NAMEN, FARBEN): if abs(zeile[1]) < 1e-9: # senkrechte Gerade x1 = c achse.axvline(grenze / zeile[0], color=farbe, linewidth=1.3, linestyle="--", alpha=0.8, label=name) else: achse.plot(x, (grenze - zeile[0] * x) / zeile[1], color=farbe, linewidth=1.3, linestyle="--", alpha=0.8, label=name) achse.set_xlim(-0.9, 11.4) achse.set_ylim(-0.9, 11) achse.set_xlabel("$x_1$ — Bots vom Typ A") achse.set_ylabel("$x_2$ — Bots vom Typ B") achse.grid(linestyle=":", alpha=0.45) for seite in ("top", "right"): achse.spines[seite].set_visible(False) def schreibe(figur, name: str) -> None: os.makedirs(BILDER, exist_ok=True) for endung in ("svg",): pfad = os.path.join(BILDER, f"{name}.{endung}") figur.savefig(pfad, format=endung, dpi=160, metadata={"Date": None} if endung == "svg" else None) print(f"geschrieben: {pfad}") plt.close(figur) def zeichne_polyeder(ecken: np.ndarray, bestes: np.ndarray) -> None: """Fundament: alle Ecken mit ihrem Zielwert - der Fundamentalsatz zum Anfassen.""" figur, achse = plt.subplots(figsize=(7.2, 5.2)) zeichne_rahmen(achse, ecken) for ecke in ecken: wert = float(C @ ecke) optimal = np.allclose(ecke, bestes) achse.scatter(*ecke, s=120 if optimal else 55, zorder=4, color=GRUEN if optimal else INDIGO) beschriftung = f"{punkt_text(ecke)}\n{euro(wert)}" if optimal: beschriftung += "\n← Optimum" weg = versatz(ecke) achse.annotate(beschriftung, ecke, textcoords="offset points", xytext=weg, fontsize=8.5, ha="left" if weg[0] >= 0 else "right", va="bottom" if weg[1] >= 0 else "top", color=GRUEN if optimal else "#334155", fontweight="bold" if optimal else "normal") achse.set_title("Der zulässige Bereich und seine Ecken\n" "$\\max\\ 150x_1 + 250x_2$", fontsize=11) achse.legend(fontsize=8, loc="upper right", frameon=False) figur.tight_layout() schreibe(figur, "kap_fundament_polyeder") def zeichne_niveaulinien(achse, besucht: list[np.ndarray]) -> None: """Eine Niveaulinie der Zielfunktion durch JEDE besuchte Ecke. Die Niveaus werden aus dem Weg gerechnet, nicht gewaehlt: Jede Ecke des Simplexwegs liegt damit auf genau einer Linie, und der Abstand der Linien ist der Fortschritt eines Schritts. Damit zeigt das Bild nicht nur, WELCHEN Weg der Simplex nimmt, sondern warum - er schiebt die Linie so weit nach aussen, wie das Polyeder es zulaesst.""" x = np.linspace(-0.9, 11.4, 100) for punkt in besucht: niveau = float(C @ punkt) achse.plot(x, (niveau - C[0] * x) / C[1], linestyle=":", color=VIOLETT, linewidth=1.3, alpha=0.9, zorder=2) # Beschriftung dort, wo die Linie im Bild noch Platz hat. Die # Nulllinie laeuft nahe am unteren Rand - ihre Beschriftung gehoert # darueber, sonst wird sie abgeschnitten. unten = niveau < 500 x_text = 1.7 if unten else 10.1 y_linie = (niveau - C[0] * x_text) / C[1] achse.text(x_text, y_linie + (0.16 if unten else -0.28), f"Z = {euro(niveau)}", fontsize=7.5, color=VIOLETT, ha="center", va="bottom" if unten else "top", zorder=6) # Der Gradient: die Richtung, in die sich die Linien verschieben lassen. start = np.array([2.6, 2.2]) ende = start + C / np.linalg.norm(C) * 2.6 achse.annotate("", xy=tuple(ende), xytext=tuple(start), arrowprops=dict(arrowstyle="-|>", color=VIOLETT, linewidth=2.2), zorder=6) achse.text(ende[0] + 0.15, ende[1], "$\\mathbf{c} = (150; 250)$", fontsize=8.5, color=VIOLETT, va="center", zorder=6) def zeichne_simplexweg(ecken: np.ndarray, bestes: np.ndarray) -> None: """LP-Kapitel: derselbe Bereich, aber der Weg der Handrechnung.""" figur, achse = plt.subplots(figsize=(7.2, 5.2)) zeichne_rahmen(achse, ecken) besucht = [np.array(p) for p in SIMPLEXWEG] zeichne_niveaulinien(achse, besucht) unbesucht = [e for e in ecken if not any(np.allclose(e, p) for p in besucht)] for ecke in unbesucht: achse.scatter(*ecke, s=55, color="#94a3b8", zorder=3) weg = versatz(ecke, mit_pfad=True) achse.annotate(punkt_text(ecke), ecke, textcoords="offset points", xytext=weg, fontsize=8, color="#64748b", ha="left" if weg[0] >= 0 else "right", va="bottom" if weg[1] >= 0 else "top") weg = np.array(besucht) achse.plot(weg[:, 0], weg[:, 1], color=AMBER, linewidth=2.4, zorder=4) for schritt, punkt in enumerate(besucht): wert = float(C @ punkt) letzter = schritt == len(besucht) - 1 achse.scatter(*punkt, s=140 if letzter else 90, zorder=5, color=GRUEN if letzter else AMBER) weg = versatz(punkt, mit_pfad=True) achse.annotate(f"{schritt}. {punkt_text(punkt)}\nZ = {euro(wert)}", punkt, textcoords="offset points", xytext=weg, fontsize=8.5, fontweight="bold", ha="left" if weg[0] >= 0 else "right", va="bottom" if weg[1] >= 0 else "top", color=GRUEN if letzter else "#334155") achse.set_title(f"Der Simplex besucht {len(besucht)} von {len(ecken)} Ecken\n" "und keine davon zweimal", fontsize=11) achse.legend(fontsize=8, loc="upper right", frameon=False) figur.tight_layout() schreibe(figur, "kap_lp_simplex_ecken") if __name__ == "__main__": ecken = berechne_ecken() werte = ecken @ C bestes = ecken[int(np.argmax(werte))] print(f"{len(ecken)} Ecken gefunden:") for ecke, wert in sorted(zip(ecken.tolist(), werte.tolist()), key=lambda p: -p[1]): print(f" ({ecke[0]:6.2f}; {ecke[1]:5.2f}) Z = {wert:9,.2f} EUR") print(f"Optimum: ({bestes[0]:g}; {bestes[1]:g})") for punkt in SIMPLEXWEG: if not any(np.allclose(np.array(punkt), e) for e in ecken): raise SystemExit(f"FEHLER: {punkt} ist keine Ecke des Polyeders.") zeichne_polyeder(ecken, bestes) zeichne_simplexweg(ecken, bestes)