{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Kapitel 11: Quadratische und nichtlineare Optimierung — KKT, Lagrange, Konvexität\n", "\n", "Begleitnotebook zu *Optimierte Entscheidungsfindung mit Python*. Die Codezellen sind identisch mit den im Buch abgedruckten Programmen.\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Einmalig ausfuehren: installiert alle im Buch verwendeten Pakete.\n", "# Lokal in einer virtuellen Umgebung genauso gueltig wie in Google Colab.\n", "%pip install --quiet ortools highspy cvxpy scipy numpy pandas polars \\\n", " scikit-learn matplotlib plotly pyomo linopy pymoo pydantic openpyxl" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Konvexität: die präzise Aussage\n", "\n", "`QP_Grundlagen.py`\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "#!/usr/bin/env python3\n", "\n", "# QP_Grundlagen.py\n", "\"\"\"\n", "Kapitel QP/NLP: Die drei Faelle aus der Konvexitaets-Tabelle (Abschnitt\n", "'Das quadratische Programm') an\n", "einem Mini-QP demonstriert: P positiv definit, P (singulaer) semidefinit,\n", "P mit negativem Eigenwert.\n", "\"\"\"\n", "\n", "import numpy as np\n", "import cvxpy as cp\n", "\n", "\n", "def loese_qp(P, q, name):\n", " n = len(q)\n", " w = cp.Variable(n)\n", " ziel = cp.Minimize(0.5 * cp.quad_form(w, P) + q @ w)\n", " bedingungen = [cp.sum(w) == 1, w >= 0]\n", " problem = cp.Problem(ziel, bedingungen)\n", "\n", " print(f\"\\n--- {name} ---\")\n", " eigenwerte = np.linalg.eigvalsh(P)\n", " print(f\"Eigenwerte von P: {np.round(eigenwerte, 4)}\")\n", " print(f\"DCP-konvex (CVXPY-Pruefung)? {problem.is_dcp()}\")\n", "\n", " if not problem.is_dcp():\n", " print(\"-> CVXPY lehnt das Problem ab, BEVOR ueberhaupt ein Solver laeuft.\")\n", " return\n", "\n", " problem.solve()\n", " print(f\"Status: {problem.status}\")\n", " print(f\"w* = {np.round(w.value, 4)}\")\n", " print(f\"Zielwert = {problem.value:.6f}\")\n", "\n", "\n", "if __name__ == \"__main__\":\n", " q = np.zeros(2)\n", "\n", " # Fall 1: P positiv definit -> eindeutiges Minimum\n", " P_definit = np.array([[2.0, 0.5], [0.5, 1.0]])\n", " loese_qp(P_definit, q, \"P positiv definit\")\n", "\n", " # Fall 2: P singulaer/semidefinit (zwei \"identische\" Assets) -> unendlich viele Minima\n", " P_semidefinit = np.array([[1.0, 1.0], [1.0, 1.0]])\n", " loese_qp(P_semidefinit, q, \"P positiv semidefinit (singulaer)\")\n", "\n", " # Fall 3: P mit negativem Eigenwert -> nicht konvex\n", " P_indefinit = np.array([[1.0, 2.0], [2.0, 1.0]])\n", " loese_qp(P_indefinit, q, \"P indefinit (negativer Eigenwert)\")\n", "\n", " print(\"\\n--- Nachweis: 'unendlich viele Minima' im semidefiniten Fall ---\")\n", " for punkt in [np.array([1.0, 0.0]), np.array([0.0, 1.0]), np.array([0.3, 0.7])]:\n", " wert = 0.5 * punkt @ P_semidefinit @ punkt\n", " print(f\" w = {punkt} -> Zielwert = {wert:.4f} (identisch, obwohl w verschieden)\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Die vier KKT-Bedingungen\n", "\n", "`KKT_Nachweis.py`\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "#!/usr/bin/env python3\n", "\n", "# KKT_Nachweis.py\n", "\"\"\"\n", "Kapitel QP/NLP: Die KKT-Bedingungen numerisch nachpruefen.\n", "\n", "Loest ein QP mit CVXPY, liest die Dualwerte aus und prueft alle vier\n", "KKT-Bedingungen einzeln nach. Das ist zugleich eine Vorlage fuer die\n", "Qualitaetssicherung eigener Modelle.\n", "\"\"\"\n", "\n", "import numpy as np\n", "import cvxpy as cp\n", "\n", "\n", "def baue_gueltige_kovarianz(vola, korrelationen, seed=0):\n", " \"\"\"\n", " Baut Sigma = D * C * D aus Volatilitaeten und einer Korrelationsmatrix.\n", " Dieses Vorgehen ist konstruktionsbedingt positiv semidefinit - im\n", " Gegensatz zum nachtraeglichen Ueberschreiben der Diagonalen, das die\n", " positive Semidefinitheit zerstoeren kann.\n", " \"\"\"\n", " C = np.array(korrelationen, dtype=float)\n", " assert np.allclose(C, C.T), \"Korrelationsmatrix muss symmetrisch sein.\"\n", " assert np.allclose(np.diag(C), 1.0), \"Diagonale der Korrelationsmatrix muss 1 sein.\"\n", " eigen = np.linalg.eigvalsh(C)\n", " assert eigen.min() > -1e-10, (\n", " f\"Korrelationsmatrix ist nicht positiv semidefinit \"\n", " f\"(kleinster Eigenwert {eigen.min():.4f}). Solche Korrelationen sind unmoeglich.\")\n", " D = np.diag(vola)\n", " return D @ C @ D\n", "\n", "\n", "if __name__ == \"__main__\":\n", " # --- Ein kleines Portfolio-QP ----------------------------------------\n", " vola = np.array([0.20, 0.14, 0.30]) # Volatilitaeten\n", " korr = [[1.00, 0.30, 0.10],\n", " [0.30, 1.00, 0.25],\n", " [0.10, 0.25, 1.00]]\n", " Sigma = baue_gueltige_kovarianz(vola, korr)\n", " mu = np.array([0.09, 0.05, 0.13]) # erwartete Renditen\n", " lam = 3.0 # Risikoaversion\n", "\n", " eigenwerte = np.linalg.eigvalsh(Sigma)\n", " print(\"=\" * 74)\n", " print(\" KKT-BEDINGUNGEN AM PORTFOLIO-QP\")\n", " print(\"=\" * 74)\n", " print(f\"Eigenwerte von Sigma: {np.round(eigenwerte, 6)}\")\n", " print(f\" -> positiv definit: {bool(eigenwerte.min() > 0)} \"\n", " f\"(Problem ist streng konvex, Loesung eindeutig)\")\n", " print(f\" -> Konditionszahl: {eigenwerte.max() / eigenwerte.min():.2f}\")\n", "\n", " # --- Modell: min lam/2 * w'Sigma w - mu'w u.d.N. sum(w)=1, w>=0 ----\n", " n = len(mu)\n", " w = cp.Variable(n)\n", " ziel = cp.Minimize(0.5 * lam * cp.quad_form(w, Sigma) - mu @ w)\n", " budget = cp.sum(w) == 1\n", " nichtnegativ = w >= 0\n", " problem = cp.Problem(ziel, [budget, nichtnegativ])\n", " problem.solve()\n", "\n", " w_opt = w.value\n", " nu = budget.dual_value # Multiplikator der Gleichung\n", " lam_i = nichtnegativ.dual_value # Multiplikatoren der Ungleichungen\n", "\n", " print(f\"\\nStatus: {problem.status}\")\n", " print(f\"Optimale Gewichte: {np.round(w_opt, 6)}\")\n", " print(f\"Zielwert: {problem.value:.6f}\")\n", " print(f\"Multiplikator der Budgetgleichung (nu): {nu:.6f}\")\n", " print(f\"Multiplikatoren der w>=0-Bedingungen: {np.round(lam_i, 6)}\")\n", "\n", " # --- KKT-Bedingungen einzeln pruefen ---------------------------------\n", " print(\"\\n--- Pruefung der vier KKT-Bedingungen ---\")\n", "\n", " # 1. Stationaritaet: grad f - lambda + nu*1 = 0\n", " # f(w) = lam/2 w'Sigma w - mu'w -> grad f = lam*Sigma w - mu\n", " # g_i(w) = -w_i <= 0 -> grad g_i = -e_i\n", " # h(w) = sum(w) - 1 = 0 -> grad h = 1\n", " grad_f = lam * (Sigma @ w_opt) - mu\n", " stationaritaet = grad_f - lam_i + nu * np.ones(n)\n", " print(f\"1. Stationaritaet : max|Residuum| = {np.abs(stationaritaet).max():.2e}\")\n", "\n", " # 2. Primale Zulaessigkeit\n", " print(f\"2. Primal zulaessig : sum(w)-1 = {w_opt.sum()-1:.2e}, \"\n", " f\"min(w) = {w_opt.min():.2e}\")\n", "\n", " # 3. Duale Zulaessigkeit\n", " print(f\"3. Dual zulaessig : min(lambda) = {lam_i.min():.2e} (muss >= 0 sein)\")\n", "\n", " # 4. Komplementaerer Schlupf: lambda_i * w_i = 0\n", " print(f\"4. Kompl. Schlupf : max|lambda_i * w_i| = \"\n", " f\"{np.abs(lam_i * w_opt).max():.2e}\")\n", "\n", " alle_ok = (np.abs(stationaritaet).max() < 1e-6\n", " and abs(w_opt.sum() - 1) < 1e-8\n", " and w_opt.min() > -1e-8\n", " and lam_i.min() > -1e-8\n", " and np.abs(lam_i * w_opt).max() < 1e-6)\n", " print(f\"\\nAlle vier KKT-Bedingungen erfuellt: {alle_ok}\")\n", "\n", " # --- Interpretation von nu -------------------------------------------\n", " print(\"\\n--- Was bedeutet nu? ---\")\n", " print(\"nu ist der Schattenpreis des Budgets: Um so viel aendert sich der\")\n", " print(\"Zielwert, wenn man statt 100 % nur 99 % investieren duerfte.\")\n", " problem2 = cp.Problem(cp.Minimize(0.5 * lam * cp.quad_form(w, Sigma) - mu @ w),\n", " [cp.sum(w) == 1.01, w >= 0])\n", " problem2.solve()\n", " print(f\" Vorhergesagt (nu * 0.01): {nu * 0.01:+.6f}\")\n", " print(f\" Tatsaechlich gemessen: {problem2.value - problem.value:+.6f}\")\n", " print(\"=\" * 74)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Praxisfall: Entropie-maximierte Kapitalallokation\n", "\n", "`Entropie_Maximierte_Allokation.py`\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "#!/usr/bin/env python3\n", "\n", "# Entropie_Maximierte_Allokation.py\n", "\"\"\"\n", "Kapitel QP/NLP: Nichtlineare Optimierung (NLP) mit scipy SLSQP.\n", "Modell: Mean-Variance-Portfolio mit Shannon-Entropie-Diversifikation.\n", "\n", "Achtung, haeufiger Fehler: Sigma per A@A.T zu erzeugen und anschliessend die\n", "Diagonale zu ueberschreiben zerstoert die positive Semidefinitheit - die\n", "Matrix kann dadurch einen negativen Eigenwert bekommen und unmoegliche\n", "Korrelationen (> 1) implizieren. Dieses Programm baut Sigma stattdessen als\n", "D * C * D aus Volatilitaeten und einer echten Korrelationsmatrix -\n", "konstruktionsbedingt immer PSD.\n", "\n", "Zusaetzlich: Multistart, weil SLSQP nur lokale Optima findet.\n", "\"\"\"\n", "\n", "import numpy as np\n", "import pandas as pd\n", "from scipy.optimize import minimize\n", "\n", "ASSETS = [\"Tech-Aktien\", \"Rohstoffe\", \"US-Treasuries\", \"Krypto\"]\n", "MU = np.array([0.14, 0.07, 0.03, 0.22]) # erwartete Jahresrenditen\n", "VOLA = np.array([0.20, 0.14, 0.07, 0.40]) # Volatilitaeten p.a.\n", "KORRELATION = np.array([\n", " [1.00, 0.25, -0.10, 0.55],\n", " [0.25, 1.00, 0.05, 0.20],\n", " [-0.10, 0.05, 1.00, -0.15],\n", " [0.55, 0.20, -0.15, 1.00],\n", "])\n", "\n", "ALPHA = 1.0 # Gewicht des Ertrags\n", "BETA = 0.015 # Gewicht der Entropie-Diversifikation\n", "UNTERGRENZE = 0.001\n", "OBERGRENZE = 0.80\n", "\n", "\n", "def baue_kovarianz() -> np.ndarray:\n", " \"\"\"Sigma = D * C * D. Immer PSD, wenn C eine gueltige Korrelationsmatrix ist.\"\"\"\n", " eig_c = np.linalg.eigvalsh(KORRELATION)\n", " if eig_c.min() < -1e-10:\n", " raise ValueError(f\"Korrelationsmatrix unmoeglich (Eigenwert {eig_c.min():.4f}).\")\n", " D = np.diag(VOLA)\n", " Sigma = D @ KORRELATION @ D\n", " eig_s = np.linalg.eigvalsh(Sigma)\n", " assert eig_s.min() > 0, \"Sigma nicht positiv definit!\"\n", " return Sigma\n", "\n", "\n", "SIGMA = baue_kovarianz()\n", "\n", "\n", "def zielfunktion(w, alpha=ALPHA, beta=BETA):\n", " \"\"\"Risiko - Ertrag - Entropiepraemie (wird minimiert).\"\"\"\n", " varianz = 0.5 * w @ SIGMA @ w\n", " ertrag = MU @ w\n", " entropie = -np.sum(w * np.log(w + 1e-12))\n", " return varianz - alpha * ertrag - beta * entropie\n", "\n", "\n", "def gradient(w, alpha=ALPHA, beta=BETA):\n", " \"\"\"Exakter analytischer Gradient - beschleunigt und stabilisiert SLSQP.\"\"\"\n", " grad_varianz = SIGMA @ w # d/dw von 0.5 w'Sigma w\n", " grad_ertrag = -alpha * MU\n", " grad_entropie = beta * (np.log(w + 1e-12) + 1.0)\n", " return grad_varianz + grad_ertrag + grad_entropie\n", "\n", "\n", "def optimiere(startpunkt, alpha=ALPHA, beta=BETA):\n", " \"\"\"alpha und beta werden durchgereicht - so bleibt die Funktion seiteneffektfrei.\"\"\"\n", " return minimize(\n", " zielfunktion, startpunkt, args=(alpha, beta), jac=gradient,\n", " method=\"SLSQP\", bounds=[(UNTERGRENZE, OBERGRENZE)] * len(MU),\n", " constraints=({\"type\": \"eq\",\n", " \"fun\": lambda w: np.sum(w) - 1.0,\n", " \"jac\": lambda w: np.ones(len(MU))}),\n", " options={\"ftol\": 1e-12, \"maxiter\": 300})\n", "\n", "\n", "if __name__ == \"__main__\":\n", " n = len(MU)\n", " eig = np.linalg.eigvalsh(SIGMA)\n", "\n", " print(\"=\" * 74)\n", " print(\" NICHTLINEARE ENTROPIE-OPTIMIERTE ASSET-ALLOKATION\")\n", " print(\"=\" * 74)\n", " print(\"Pruefung der Kovarianzmatrix:\")\n", " print(f\" Eigenwerte: {np.round(eig, 6)}\")\n", " print(f\" Positiv definit: {bool(eig.min() > 0)}\")\n", " print(f\" Groesste Korrelation ausserhalb der Diagonale: \"\n", " f\"{np.abs(KORRELATION - np.eye(n)).max():.2f} (muss <= 1 sein)\")\n", "\n", " # --- Multistart: SLSQP findet nur lokale Optima -----------------------\n", " rng = np.random.default_rng(7)\n", " startpunkte = [np.ones(n) / n] # Gleichgewichtung\n", " for _ in range(9):\n", " z = rng.random(n) + 0.05\n", " startpunkte.append(z / z.sum())\n", "\n", " ergebnisse = [optimiere(s) for s in startpunkte]\n", " erfolgreich = [r for r in ergebnisse if r.success]\n", " bestes = min(erfolgreich, key=lambda r: r.fun)\n", " zielwerte = np.array([r.fun for r in erfolgreich])\n", "\n", " print(f\"\\nMultistart mit {len(startpunkte)} Startpunkten:\")\n", " print(f\" Erfolgreich konvergiert: {len(erfolgreich)}\")\n", " print(f\" Spannweite der Zielwerte: {zielwerte.max() - zielwerte.min():.2e}\")\n", " print(\" -> \" + (\"alle Startpunkte fuehren zum selben Optimum (Indiz fuer Konvexitaet)\"\n", " if zielwerte.max() - zielwerte.min() < 1e-6\n", " else \"ACHTUNG: verschiedene lokale Optima gefunden!\"))\n", "\n", " w_opt = bestes.x\n", " rendite = MU @ w_opt\n", " volatilitaet = np.sqrt(w_opt @ SIGMA @ w_opt)\n", " entropie = -np.sum(w_opt * np.log(w_opt))\n", "\n", " print(f\"\\nKonvergenz: {bestes.message}\")\n", " print(f\"Iterationen: {bestes.nit}\")\n", " print(f\"Erwartete Jahresrendite: {rendite * 100:6.2f} %\")\n", " print(f\"Erwartete Volatilitaet: {volatilitaet * 100:6.2f} %\")\n", " print(f\"Diversifikations-Entropie: {entropie:.4f} \"\n", " f\"(Maximum bei Gleichgewichtung: {np.log(n):.4f})\")\n", "\n", " print(\"\\n\" + pd.DataFrame({\n", " \"Asset\": ASSETS,\n", " \"Erw. Rendite\": [f\"{r*100:.1f} %\" for r in MU],\n", " \"Volatilitaet\": [f\"{v*100:.1f} %\" for v in VOLA],\n", " \"Gewicht\": [f\"{w*100:6.2f} %\" for w in w_opt],\n", " }).to_string(index=False))\n", "\n", " # --- Vergleich: was passiert ohne Entropieterm? ----------------------\n", " ohne = optimiere(np.ones(n) / n, beta=0.0)\n", " print(\"\\n--- Wirkung des Entropieterms ---\")\n", " print(f\" Mit Entropie (beta={BETA}): Gewichte {np.round(w_opt * 100, 1)}\")\n", " print(f\" Ohne Entropie (beta=0): Gewichte {np.round(ohne.x * 100, 1)}\")\n", " print(\" Der Entropieterm zieht Kapital aus der Spitzenposition heraus,\")\n", " print(\" ohne dass eine harte Obergrenze noetig waere.\")\n", " print(\"=\" * 74)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Was ein lokales Optimum praktisch bedeutet\n", "\n", "`Lokale_Optima_Multistart.py`\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "#!/usr/bin/env python3\n", "\n", "# Lokale_Optima_Multistart.py\n", "\"\"\"\n", "Kapitel QP/NLP: Was passiert, wenn die Konvexitaet fehlt.\n", "\n", "CVXPY verweigert nicht-konvexe Probleme - das ist sein Schutzmechanismus.\n", "scipy.optimize.minimize verweigert nichts. Es rechnet, meldet 'success: True'\n", "und liefert ein Ergebnis. Nur ist das Ergebnis dann kein Optimum, sondern\n", "irgendein lokales Minimum, das vom Startpunkt abhaengt.\n", "\n", "Beispiel aus dem Einkauf: 700 Tonnen Rohstoff werden auf vier Lieferanten\n", "verteilt. Jeder gewaehrt einen MENGENRABATT - der Stueckpreis faellt, je mehr\n", "man bei ihm bestellt:\n", "\n", " preis_i(q) = basis_i * (1 - rabatt_i * (1 - exp(-q / skala_i)))\n", "\n", "Genau das macht die Zielfunktion nicht-konvex: Grosse Bestellungen lohnen\n", "sich ueberproportional, es gibt also mehrere sinnvolle \"Cluster\"-Loesungen -\n", "und dazwischen schlechtere Taeler.\n", "\n", "Das Programm zeigt drei Dinge:\n", " 1. Ein einzelner Lauf liefert ein plausibles Ergebnis - ohne jede Warnung.\n", " 2. 200 Startpunkte foerdern mehrere verschiedene lokale Optima zutage.\n", " 3. Der Unterschied zwischen bestem und schlechtestem betraegt hier 8 %.\n", "\n", "Benoetigt: numpy, scipy\n", "\"\"\"\n", "\n", "from __future__ import annotations\n", "\n", "import numpy as np\n", "from scipy.optimize import minimize\n", "\n", "# Vier Lieferanten\n", "NAMEN = [\"Nord AG\", \"Ost GmbH\", \"Sued KG\", \"West SE\"]\n", "BASISPREIS = np.array([50.0, 47.0, 53.0, 45.0]) # EUR je Tonne ohne Rabatt\n", "MAX_RABATT = np.array([0.30, 0.18, 0.35, 0.12]) # hoechstmoeglicher Rabatt\n", "RABATT_SKALA = np.array([120.0, 260.0, 90.0, 300.0]) # wie schnell er greift\n", "KAPAZITAET = np.array([400.0, 400.0, 400.0, 400.0])\n", "BEDARF = 700.0\n", "\n", "RNG = np.random.default_rng(0)\n", "\n", "\n", "def gesamtkosten(menge: np.ndarray) -> float:\n", " \"\"\"Einkaufskosten bei mengenabhaengigem Stueckpreis.\n", "\n", " Der Rabatt waechst mit der Bestellmenge und laeuft gegen MAX_RABATT.\n", " Dadurch ist der Stueckpreis fallend - und die Gesamtkostenfunktion\n", " nicht mehr konvex.\n", " \"\"\"\n", " stueckpreis = BASISPREIS * (1 - MAX_RABATT * (1 - np.exp(-menge / RABATT_SKALA)))\n", " return float(stueckpreis @ menge)\n", "\n", "\n", "NEBENBEDINGUNGEN = [{\"type\": \"eq\", \"fun\": lambda q: q.sum() - BEDARF}]\n", "GRENZEN = [(0.0, k) for k in KAPAZITAET]\n", "\n", "\n", "def optimiere_von(startpunkt: np.ndarray):\n", " \"\"\"Ein Lauf von einem gegebenen Startpunkt aus.\"\"\"\n", " return minimize(gesamtkosten, startpunkt, method=\"SLSQP\",\n", " bounds=GRENZEN, constraints=NEBENBEDINGUNGEN)\n", "\n", "\n", "def zufaelliger_start() -> np.ndarray:\n", " \"\"\"Zufaellige Aufteilung, die den Bedarf bereits erfuellt.\"\"\"\n", " anteil = RNG.random(len(NAMEN))\n", " return anteil / anteil.sum() * BEDARF\n", "\n", "\n", "def zeige_plan(titel: str, menge: np.ndarray, kosten: float) -> None:\n", " print(f\"\\n{titel}\")\n", " for name, m in zip(NAMEN, menge):\n", " anteil = m / BEDARF * 100\n", " print(f\" {name:<10} {m:7.1f} t ({anteil:4.1f} %)\")\n", " print(f\" {'Gesamtkosten':<10} {kosten:9,.2f} EUR\")\n", "\n", "\n", "if __name__ == \"__main__\":\n", " print(\"=\" * 78)\n", " print(\" NICHT-KONVEX: DERSELBE CODE, VERSCHIEDENE ERGEBNISSE\")\n", " print(\"=\" * 78)\n", " print(f\"{BEDARF:.0f} t Rohstoff auf {len(NAMEN)} Lieferanten mit Mengenrabatt.\")\n", "\n", " # --- 1. Ein einziger Lauf, so wie man es zuerst schreibt --------------\n", " erster = optimiere_von(np.full(len(NAMEN), BEDARF / len(NAMEN)))\n", " print(f\"\\n[1] EIN Lauf, Startpunkt 'gleichmaessig verteilt'\")\n", " print(f\" scipy meldet: success={erster.success}, \"\n", " f\"'{erster.message}'\")\n", " zeige_plan(\" Ergebnis:\", erster.x, erster.fun)\n", " print(\"\\n Nichts an dieser Ausgabe deutet darauf hin, dass etwas fehlt.\")\n", "\n", " # --- 1b. Der kaufmaennisch naheliegende Startpunkt --------------------\n", " # \"Kaufe bei den beiden Lieferanten mit dem guenstigsten Basispreis\" -\n", " # West SE (45) und Ost GmbH (47), jeweils bis zur Kapazitaetsgrenze.\n", " guenstigste = np.argsort(BASISPREIS)[:2]\n", " kaufmaennisch = np.zeros(len(NAMEN))\n", " rest = BEDARF\n", " for i in guenstigste:\n", " kaufmaennisch[i] = min(KAPAZITAET[i], rest)\n", " rest -= kaufmaennisch[i]\n", " zweiter = optimiere_von(kaufmaennisch)\n", " print(f\"\\n[1b] EIN Lauf, Startpunkt 'die zwei mit dem guenstigsten Basispreis'\")\n", " print(f\" scipy meldet: success={zweiter.success}\")\n", " zeige_plan(\" Ergebnis:\", zweiter.x, zweiter.fun)\n", " print(f\"\\n Dasselbe Programm, derselbe Aufruf, ein anderer Startpunkt -\")\n", " print(f\" und {zweiter.fun - erster.fun:,.2f} EUR Unterschied \"\n", " f\"({(zweiter.fun / erster.fun - 1) * 100:.1f} %).\")\n", "\n", " # --- 2. Multistart: dasselbe Problem, viele Startpunkte ---------------\n", " laeufe = []\n", " for _ in range(200):\n", " ergebnis = optimiere_von(zufaelliger_start())\n", " if ergebnis.success:\n", " laeufe.append((float(ergebnis.fun), ergebnis.x))\n", "\n", " # Ergebnisse, die sich um weniger als 1 Cent unterscheiden, sind dasselbe\n", " # lokale Optimum - zusammenfassen, sonst zaehlt man Rundungsrauschen.\n", " optima: list[tuple[float, np.ndarray]] = []\n", " for wert, plan in sorted(laeufe, key=lambda t: t[0]):\n", " if not optima or abs(wert - optima[-1][0]) > 0.01:\n", " optima.append((wert, plan))\n", "\n", " print(\"\\n\" + \"=\" * 78)\n", " print(f\"[2] 200 zufaellige Startpunkte -> {len(laeufe)} erfolgreiche Laeufe\")\n", " print(f\" darunter {len(optima)} VERSCHIEDENE lokale Optima:\")\n", " print()\n", " print(f\" {'Rang':>5} {'Kosten':>13} {'Abstand zum besten':>20} Aufteilung (t)\")\n", " print(\" \" + \"-\" * 70)\n", " bester = optima[0][0]\n", " for rang, (wert, plan) in enumerate(optima, start=1):\n", " abstand = (wert / bester - 1) * 100\n", " aufteilung = \" \".join(f\"{m:5.0f}\" for m in plan)\n", " print(f\" {rang:>5} {wert:>13,.2f} {abstand:>19.2f} % {aufteilung}\")\n", "\n", " # --- 3. Was das kostet ------------------------------------------------\n", " schlechtester = optima[-1]\n", " print(\"\\n\" + \"=\" * 78)\n", " print(\" WAS AUF DEM SPIEL STEHT\")\n", " print(\"=\" * 78)\n", " zeige_plan(\"Bester gefundener Plan:\", optima[0][1], optima[0][0])\n", " zeige_plan(\"Schlechtestes lokales Optimum:\", schlechtester[1], schlechtester[0])\n", " unterschied = schlechtester[0] - optima[0][0]\n", " print(f\"\\n Unterschied: {unterschied:,.2f} EUR \"\n", " f\"({unterschied / optima[0][0] * 100:.1f} %)\")\n", " print(f\"\\n Lauf [1] (gleichmaessiger Start): {erster.fun:>10,.2f} EUR\")\n", " print(f\" Lauf [1b] (kaufmaennischer Start): {zweiter.fun:>10,.2f} EUR\"\n", " f\" <- {(zweiter.fun / optima[0][0] - 1) * 100:.1f} % ueber dem besten\")\n", " print()\n", " print(\" Bemerkenswert: Der kaufmaennisch NAHELIEGENDE Startpunkt fuehrt in\")\n", " print(\" das schlechteste Ergebnis von allen - schlechter als jedes der 200\")\n", " print(\" zufaellig gefundenen lokalen Optima. Wer beim guenstigsten\")\n", " print(\" Basispreis anfaengt, uebersieht, dass hier der Mengenrabatt\")\n", " print(\" entscheidet und nicht der Listenpreis.\")\n", "\n", " print()\n", " print(\"Drei Konsequenzen fuer die Praxis:\")\n", " print(\" 1. 'success: True' heisst bei nicht-konvexen Problemen NICHT 'optimal'.\")\n", " print(\" Es heisst nur: 'Ich bin an einer Stelle angekommen, an der es in\")\n", " print(\" keine Richtung mehr bergab geht.'\")\n", " print(\" 2. Ein einzelner Lauf ist wertlos. Nehmen Sie viele Startpunkte und\")\n", " print(\" berichten Sie die STREUUNG mit - sie ist Ihre einzige Auskunft\")\n", " print(\" darueber, wie zerklueftet die Landschaft ist.\")\n", " print(\" 3. Auch Multistart liefert KEINE Garantie. Dass hier nichts unter\")\n", " print(f\" {bester:,.2f} EUR gefunden wurde, beweist nicht, dass es nichts gibt.\")\n", " print(\"=\" * 78)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "name": "python", "version": "3.11" } }, "nbformat": 4, "nbformat_minor": 5 }