operations_research/bilder_04/erzeuge_entartung_polyeder.py

#!/usr/bin/env python3

# erzeuge_entartung_polyeder.py
"""
Erzeugt das Entartungs-Diagramm zum LP-Kapitel:

    bilder_04/kap_lp_entartung.svg

Links die Geometrie: Drei Geraden laufen durch DIESELBE Ecke. In zwei
Dimensionen legen schon zwei eine Ecke fest - eine ist zu viel, und genau das
ist Entartung. Rechts die Folge davon: Der Schattenpreis jeder Nebenbedingung
ist nicht mehr eine Zahl, sondern eine Spanne.

Warum beides in ein Bild gehoert: Die Geometrie allein erklaert nicht, warum
zwei Solver verschiedene Dualwerte melden duerfen und beide recht haben. Die
Spannen allein sehen nach einem Fehler aus. Nebeneinander ergibt es Sinn.

Modell und Spannen kommen aus Toleranzen_und_Entartung.py - importiert, nicht
abgeschrieben. pruefe_ecke() rechnet zusaetzlich nach, dass wirklich alle drei
Geraden durch den Punkt laufen; waere das nicht so, zeigte das Bild etwas
anderes als der Text daneben behauptet.

Aufruf (aus dem Repository-Wurzelverzeichnis):
    python3 bilder_04/erzeuge_entartung_polyeder.py

Benoetigt: numpy, scipy, matplotlib
"""

from __future__ import annotations

import os
import sys

import numpy as np

HIER = os.path.dirname(os.path.abspath(__file__))
BASIS = os.path.dirname(HIER)
sys.path.insert(0, HIER)
sys.path.insert(0, os.path.join(
    BASIS, "Operations_Research_mit_Python_Version_04_Programme"))

from stil_04 import FARBEN, speichere                        # noqa: E402
import matplotlib.pyplot as plt                              # noqa: E402
from Toleranzen_und_Entartung import (A_UB, B_UB, NAMEN,     # noqa: E402
                                      schattenpreis_spanne)

ECKE = np.array([4 / 3, 4 / 3])
LINIENFARBEN = (FARBEN["haupt"], FARBEN["zweit"], FARBEN["akzent"])


def pruefe_ecke() -> None:
    """Alle drei Nebenbedingungen muessen in der Ecke aktiv sein."""
    schlupf = B_UB - A_UB @ ECKE
    aktiv = np.abs(schlupf) < 1e-9
    if not aktiv.all():
        raise SystemExit(
            f"Entartungsbild: In {tuple(np.round(ECKE, 3))} sind nur "
            f"{int(aktiv.sum())} von 3 Bedingungen aktiv (Schlupf "
            f"{np.round(schlupf, 6)}). Ohne drei aktive Geraden gibt es keine "
            f"Entartung zu zeigen.")


def zeichne(spannen) -> None:
    figur, (links, rechts) = plt.subplots(1, 2, figsize=(9.4, 4.3),
                                          gridspec_kw={"width_ratios": [1, 1]})

    # --- links: die drei Geraden ---------------------------------------
    x = np.linspace(0, 3.2, 200)
    for (a, b), grenze, name, farbe in zip(A_UB, B_UB, NAMEN, LINIENFARBEN):
        y = (grenze - a * x) / b
        links.plot(x, y, "-", color=farbe, linewidth=2.0, label=name, zorder=3)

    # Zulaessiger Bereich: Schnitt aller Halbraeume, plus x >= 0.
    gitter = np.linspace(0, 3.2, 400)
    xx, yy = np.meshgrid(gitter, gitter)
    zulaessig = np.ones_like(xx, dtype=bool)
    for (a, b), grenze in zip(A_UB, B_UB):
        zulaessig &= (a * xx + b * yy <= grenze + 1e-12)
    links.contourf(xx, yy, zulaessig.astype(float), levels=[0.5, 1.5],
                   colors=[FARBEN["gut"]], alpha=0.13, zorder=1)

    links.plot(*ECKE, "o", color=FARBEN["fehler"], markersize=9, zorder=5)
    links.annotate(r"Optimum $(4/3,\ 4/3)$" "\n" "drei Geraden, eine Ecke",
                   xy=tuple(ECKE), xytext=(18, 26),
                   textcoords="offset points", fontsize=9,
                   color=FARBEN["fehler"], fontweight="bold",
                   arrowprops=dict(arrowstyle="->", color=FARBEN["fehler"],
                                   linewidth=1.2))

    links.set_xlim(0, 3.0)
    links.set_ylim(0, 3.0)
    links.set_xlabel(r"$x_1$")
    links.set_ylabel(r"$x_2$")
    links.set_title("Eine Bedingung zu viel", fontsize=10.5,
                    color=FARBEN["text"])
    links.legend(fontsize=8.5, loc="upper right", frameon=False)

    # --- rechts: die Spanne je Schattenpreis ---------------------------
    y_pos = np.arange(len(NAMEN))[::-1]
    for stelle, ((unten, oben), name, farbe) in enumerate(
            zip(spannen, NAMEN, LINIENFARBEN)):
        y = y_pos[stelle]
        rechts.plot([unten, oben], [y, y], "-", color=farbe, linewidth=7,
                    solid_capstyle="round", alpha=0.75, zorder=3)
        rechts.plot([unten, oben], [y, y], "|", color=farbe, markersize=14,
                    markeredgewidth=2, zorder=4)
        rechts.text(oben + 0.012, y, f"{unten:.2f} – {oben:.2f} €",
                    va="center", fontsize=9, color=FARBEN["text"])
        rechts.text(-0.012, y, name, va="center", ha="right", fontsize=9,
                    color=FARBEN["text"])

    rechts.set_yticks([])
    rechts.set_ylim(-0.7, len(NAMEN) - 0.3)
    rechts.set_xlim(-0.16, max(o for _, o in spannen) * 1.9)
    rechts.set_xlabel("Schattenpreis in € je zusätzlicher Stunde")
    rechts.set_title("Was der Solver melden darf — jeder Wert der Spanne "
                     "ist richtig", fontsize=10.5, color=FARBEN["text"])
    rechts.grid(axis="y", visible=False)
    rechts.spines["left"].set_visible(False)

    speichere(figur, "kap_lp_entartung")


if __name__ == "__main__":
    pruefe_ecke()
    zielwert = float(np.ones(2) @ ECKE)
    spannen = schattenpreis_spanne(zielwert)
    for name, (unten, oben) in zip(NAMEN, spannen):
        print(f"  {name:<18} {unten:>7.4f} bis {oben:>7.4f} EUR")
    zeichne(spannen)