Lineare Optimierung

Lineare Optimierung
Lineare Optimierung

Das Folgende setzt Mengenlehre, Topologie und Nichtstandardanalysis voraus.

Durchmessersatz für Polytope und Polyeder: Jeder Durchmesser eines durch \(m\) Restriktionen gegebenen \(n\)-dimensionalen Polytops bzw. Polyeders mit \(m, n \in {}^{\omega}\mathbb{N}_{\ge 2}\) ist maximal \(2(m + n – 3)\).

Beweis: Maximal \(\acute{m}\) Hyperebenen lassen sich zu einem unvollständigen Zyklus der Dimension 2 zusammenstellen und es gibt maximal \(n – 2\) Ausweichmöglichkeiten zur Seite in den restlichen Dimensionen. Das Überwinden jeder minimalen Strecke benötigt maximal zwei Kanten und liefert den Faktor 2.\(\square\)

Definition: Sei \(\vartheta := \omega\,{}_{\epsilon}\omega\) und \({|| \cdot ||}_{1}\) die 1-Norm. Geht \(x \in X \subseteq {}^{\omega}\mathbb{R}^{\omega}\) in den nächsten Schritt über, folgt ihm \({}^{*}\) mit \(\Delta x := {x}^{*} – x\). Der Einheitsvektor von \(x\) ist \({_1}{x} := x/||x||\), wobei \(_1{0}\) undefiniert ist. Seien \(x, y \in X \) in \((x, y, …)^T\) Zeilenvektoren. Benötigt ein Verfahren Rechenzeit in Sekunden und Speicher in Bits von \(\mathcal{O}({\omega}^{\mathcal{O}(1)})\), ist es polynomiell. Ist eine der beiden Größen \(\mathcal{O}({\epsilon}^{|\mathcal{O}(\omega)|})\), ist es exponentiell. Das Eigenprodukt (bisher: Determinante) einer quadratischen Matrix sei das Produkt ihrer Eigenwerte.\(\triangle\)

Satz: Das Simplexverfahren1Dantzig, George B.: Lineare Programmierung und Erweiterungen; 1. Aufl.; 1966; Springer; Berlin. löst fast immer ein LP in \(\mathcal{O}({\omega}^3\vartheta)\).

Beweis und Algorithmus: Sei \(P := \{x \in {}^{\omega}\mathbb{R}^{n} : Ax \le b, b \in {}^{\omega}\mathbb{R}^{m}, A \in {}^{\omega}\mathbb{R}^{m \times n}, m, n \in {}^{\omega}\mathbb{N}^{*}\}\) der zulässige Bereich des LPs max \(\{{c}^{T}x : c \in {}^{\omega}\mathbb{R}^{n}, x \in P\}\). Dessen Dualisierung erzielt \({x}^{*} \in {}^{\omega}\mathbb{R}_{\ge 0}^{m}\). Ein (großes) \(x^+ \in {}^{\omega}\mathbb{R}_{\geq0}\) erreicht \(x^* \in {}^{\omega}\mathbb{R}_{\geq 0}^n\) mit \(x := x^* – {x^+}{\upharpoonright}_n\). Gilt \(b \ge 0\) nicht, liefert das LP max \(\{-h \in {}^{\omega}\mathbb{R}_{\le 0} : Ax – b \le h{\upharpoonright}_m\}\) für den Startpunkt \(x_0 := 0\) ein \(x \in {P}_{\ge 0}\). Start- und Zielwert von \(h\) sei \(|\text{min } \{{b}_{1}, …, {b}_{m}\}|\) bzw. 0. Nur ein Pivotschritt bewirkt \({b}^{*} \ge 0\). Seien \(i, j, k \in {}^{\omega}\mathbb{N}^{*}\) und sei \({a}_{i}^{T}\) der \(i\)-te Zeilenvektor von \(A\). Gilt \({c}_{j} \le 0\) für alle \(j\), ist das LP gelöst.

Gibt es ein \({c}_{j} > 0\) mit \({a}_{ij} \le 0\) für alle \(i\), ist das LP positiv unbeschränkt. Existiert ein \({c}_{j} = 0\) mit \({a}_{ij} \le 0\) für alle \(i\), lassen sich \({c}_{j}\) und \({A}_{.j}\) wie die \({b}_{i}\) und \({a}_{i}\) mit \({a}_{ij} < 0\) vorerst entfernen. Unlösbar ist auch \({a}_{ij}{x}_{j} \ge 0 > {b}_{i}\) für alle \(j\). Das Teilen aller \({a}_{i}^{T}x \le {b}_{i}\) durch \(||{a}_{i}||\) und aller \({c}_{j}\) bzw. \({a}_{ij}\) durch min \(|{a}_{ij}|\) für \({a}_{ij} \ne 0\) sowie jedes \(j\) ist ggf. zu wiederholen und am Ende rückgängig zu machen. Restriktionen mit \({a}_{i} \le 0\) werden ggf. entfernt.

Die Wahl eines \({b_k}/{a}_{kj} :=\) min \(\{{b_i}/{a}_{ij} : {a}_{ij} > 0\}\) für jedes \({c}_{j} > 0\) sowie Nichtbasisvariablen \({x}_{j}\) liefert mit \({x}_{j}^{*} = {x}_{j} + {b}_{k}/{a}_{kj}\) ein nächstes \({x}^{*} \in P_{\ge 0}\). Bei Wahl der steilsten Kante bildet \({a}_{kj}\) das Pivotelement zu dem \({x}_{j}\), das \({c}^{T}{_1}{\Delta x}\) bzw. \({c}_{j}^{2}/\grave{v}_{j}^{2}\) mit \(v_j := ||{A}_{.j}||\) in der \(k\)-ten Restriktion maximiert. Mehrere Maxima erlauben die Pivotregel des Bestwerts max\({}_{k,j} {c}_{j}{b}_{k}/{a}_{kj}\) zu verwenden oder (langsamer) des kleinsten Winkels min \({{_{(1)}}1{\upharpoonright}^{T}_n}\;{_1}{c}^{*}\).

Wiederholt sich \(c^Tx^* = c^Tx\), werden die Restriktionen mit \({b}_{i} = 0\) um ein \(d \in {}^{\omega}\mathbb{R}_{>0}\) abgeschwächt. Mit \({b}_{i} := ||{a}_{i}||\) steht \(d\) ggf. nicht im Tableau. Liegt eine mehrfache Ecke ausnahmsweise erneut vor, wird zu jedem \({b}_{i}\) einfach \(||{a}_{i}||\) addiert. Sie zu verlassen kann das Lösen eines LPs mit \(c > 0\) und \(b = 0\) erfordern. Danach werden die Abschwächungen (Perturbationen) zurückgenommen. Diese Aufgabe kann sich auch am Ende stellen, wenn weitere Lösungen des LPs zu bestimmen und mindestens zwei \(c_j = 0\) sind.

Dies vermeidet die sich stark unterscheidenden Abschwächungen bspw. der lexikographischen Methode. Die Rechteckregel erlaubt \({c}_{j}^{*}, {a}_{ij}^{*}\) und \({b}_{i}^{*}\) einfach zu berechnen2vgl. Vanderbei, Robert J.: Linear Programming; 3rd Ed.; 2008; Springer; New York, S. 63.. Auch die Regel von Bland ist nicht optimal. Die evtl. exponentielle Laufzeit in der linearen Optimierung lässt sich nicht verhindern: Sind Polytope z. B. vom Klee-Minty- oder Jeroslow-Typ, unterscheiden sich ihre (entfernten) Ecken numerisch kaum. In Übereinstimmung mit dem Forschungsstand folgt die Behauptung.\(\square\)

Bit-Quadrat-Satz (kurz: Bit\({}^2\)-Satz): Wiederholungen erlauben mehr als \(\mathcal{O}(n^2)\) Multiplikationen und Additionen eines Algorithmus von Zahlen der Bitlänge \(m \le n\) mit \(m, n \in {}^{\omega}\mathbb{N}^{*}\) auf \(\mathcal{O}(n^2)\) reduzieren.\(\square\)

Beispiele: Das Matrizenprodukt \((c_{ij}) := \left({\LARGE{\textbf{+}}}_{k=1}^n{a_{ik}b_{kj}}\right)\) von \((a_{ij}), (b_{ij}) \in {}^{\omega}\mathbb{C}^{n \times n}\) über Produkte \({\LARGE{\textbf{+}}}_{k=1}^n{b_{jk}4^{km}}\) mit \(a_{ij}\) und \(i, j \in {}^{\omega}\mathbb{N}_{\le n}^{*}\), die Bestimmung einer \(n\)-dimensionalen Taylorreihe (s. Nichtstandardanalysis), der Gauß-Jordan-Algorithmus, das Simplex- und das Gram-Schmidtsche Orthogonalisierungsverfahren.

Bemerkung: Ein zu großes \(m\) wird rechnerisch geeignet herabgesetzt. Mit Simultanvergleichen ermöglicht Quaternionenkodierung Matrizen gemeinsam in der Cloud zu multiplizieren.

Satz: Das Intexverfahren löst die meisten lösbaren LPs mit Inter-/Extrapolationen in \(\mathcal{O}({\vartheta}^3)\).

Beweis und Algorithmus: Sei \(z := m + n\) und \(d \in [0, 1]\) die Dichte von \(A\). Zuerst werden \({b}^{T}y – {c}^{T}x \le 0, Ax \le b\) sowie \({A}^{T}y \ge c\) normiert und skaliert. \(P_r := \{(x, y)^T \in {}^{\omega}\mathbb{R}_{\ge 0}^{z} : {b}^{T}y – {c}^{T}x \le r \in [0, \rho], Ax – b \le r{\upharpoonright}_m, c – {A}^{T}y \le r{\upharpoonright}_n\}\) habe den Radius \(\rho := s|\min \; \{b_1, …, b_m, -c_1, …, -c_n\}|\) und den Skalierungsfaktor \(s \in [1, 2]\).

Es folgt \(0{\upharpoonright}_{z} \in \partial P_{\rho}\). Nach dem starken Dualitätssatz3a. a. O., S. 60 – 65. löst das LP min \(\{ r \in [0, \rho] : (x, y)^T \in P_r\}\) ebenso die LPs max \(\{{c}^{T}x : c \in {}^{\omega}\mathbb{R}^{n}, x \in {P}_{\ge 0}\}\) und min \(\{{b}^{T}y : y \in {}^{\omega}\mathbb{R}_{\ge 0}^{m}, {A}^{T}y \ge c\}\). Seine Lösung ist der geometrische Schwerpunkt \(g\) des Polytops \(P_0\). Mit \(p_k^* := \text{min}\,\check{p}_k + \text{max}\,\check{p}_k\) für \(k = 1, …, \grave{z}\) wird \(g\) durch \(p_0 := (x_0, y_0, r_0)^T\) approximiert, bis \(||\Delta p||_1\) hinreichend klein ist. Die Lösung \(t^o(x^o, y^o, r^o)^T\) des zweidimensionalen LPs min \(\{ r \in [0, \rho] : t \in {}^{\omega}\mathbb{R}_{> 0}, t(x_0, y_0)^T \in P_r\}\) approximiert \(g\) besser.

Es wird \(r \le \check{\rho}\) erreicht und mit \(t^o(x^o, y^o)^T\) wiederholt, bis ggf. \(g \in P_0\) in \(\mathcal{O}({}_2\rho^2dmn)\) berechnet ist. Das Lösen aller zweidimensionalen LPs \(\text{min}_k r_k\) durch Bisektionsverfahren für \(r_k \in {}^{\omega}\mathbb{R}_{\ge 0}\) und \(k = 1, …, z\) in jeweils \(\mathcal{O}({\vartheta}^2)\) ermittelt \(q \in {}^{\omega}\mathbb{R}^k\) mit \(q_k := \Delta p_k \Delta r_k/r\) und \(r := \text{min}_k \Delta r_k\). Vereinfacht sei \(|\Delta p_1| = … = |\Delta p_{z}|\). Hierbei wäre min \(r_{\grave{z}}\) für \(p^* := p + wq\) mit \(w \in {}^{\omega}\mathbb{R}_{\ge 0}\) ebenso zu lösen. Folgt \(\text{min}_k \Delta r_k r = 0\), wird die Berechnung gestoppt, andernfalls wiederholt bis min \(r = 0\) oder min \(r > 0\) feststeht.\(\square\)

Bemerkung: Ist kein Missbrauch für intransparente oder schlechte Entscheidungen zu befürchten, werden Implementierungsdetails veröffentlicht. Simplexverfahren und Face-Algorithmus4Pan, Ping-Qi: Linear Programming Computation; 1st Ed.; 2014; Springer; New York, S. 580 f. lösen LPs für kleine \(m\) und \(n\) schneller. Das Intexverfahren übertrifft die bekannten (Worst-Case-)LP-Algorithmen in \(\mathcal{O}(\omega^{37/18}\vartheta)\). Absichtliche Rundungsfehler in \(b\) und \(c\) machen jede Lösung eindeutig.

Bemerkung: Falls erforderlich schwächt ein gleicher kleiner Betrag die Restriktionen vorübergehend ab. Die Laufzeit für \(m, n \in {}^{\nu}\mathbb{N}^*\) ist \(\mathcal{O}(dmn)\). Wird es für verteiltes Rechnen in \({}^{\nu}\mathbb{R}^{\nu}\) optimiert, beträgt sie nur \(\mathcal{O}(1)\). Gilt \(\acute{v}_0 \le \lfloor v_j \rfloor – v_j\) für ganzzahlige Lösungen \(v_j\), bewährt es sich gut für (gemischt-) ganzzahlige Probleme und die (nicht-)konvexe (Pareto-)Optimierung (nach der Natur5Vasuki, A: Nature-Inspired Optimization Algorithms; 1st Ed.; 2020; CRC Press; Boca Raton.).

Bemerkung: Es macht sich das Zusammenfallen des berechneten geometrischen Schwerpunkts und des Inkugelmittelpunkts in \(P_0\) zunutze. Das (paarweise) Drehen der Koordinatenachsen des ersteren wie von \(c\) durch eine einfache Rotationsmatrix sichert dies in vergleichbarer Laufzeit ab. Zahlen der Länge \(\mathcal{O}({\omega})\) lassen sich bekanntlich nur in \(\mathcal{O}(\vartheta)\) abarbeiten. Auch das Folgende lässt sich einfach ins Komplexe überführen.

Folgerung: Das LP max \(\{{||x – {x}^{o}||}_1 : {c}^{T}x = {c}^{T}{x}^{o}, Ax \le b, x – {x}^{o} \in {[-1, 1]}^{n}, x \in {}^{\omega}\mathbb{R}_{\ge 0}^{n}\}\) kann zur ersten Lösung \({x}^{o}\) ggf. eine zweite in \(\mathcal{O}({\vartheta}^3)\) bestimmen, wobei \({y}^{o}\) analog behandelt werden kann.\(\square\)

Folgerung: Das LP max \(\{\nu|\lambda_j| + ||x_j||_1: Ax_j = \lambda_j x_j, |\lambda_j| \in [0, r_j], x_j \in {[-1, 1]}^{n*}\}\) kann für jedes \(j = 1, …, n\) den Eigenwert \(\lambda_j \in {}^{\omega}\mathbb{R}\) und den Eigenvektor \(x_j \in {}^{\omega}\mathbb{R}^{n}\) der Matrix \(A \in {}^{\omega}\mathbb{R}^{n \times n}\) in \(\mathcal{O}({\vartheta}^3)\) bestimmen.\(\square\)

Folgerung: Die LPs max \(\{{x}_{j} : Ax = 0\}\) liefern alle Lösungen von \(Ax = b\) in \(\mathcal{O}({\vartheta}^3)\). Die Matrix \(A\) ist genau dann regulär, wenn das LP max \(\{{||x||}_{1} : Ax = 0\} = 0\) ist.\(\square\)

Folgerung: Sei \({\alpha }_{j} := {A}_{.j}^{-1}\) für \(j = 1, …, n\) zur Matrix \({A}^{-1} \in {}^{\omega}\mathbb{R}^{n \times n}\) und sei \({\delta}_{ij}\) das Kronecker-Delta. Für reguläres \(A\) ist das Eigenprodukt \(\ne 0\) und lässt sich jedes \({A \alpha }_{j} = {({\delta}_{1j}, …, {\delta}_{nj})}^{T}\) in \(\mathcal{O}({\vartheta}^3)\) lösen.\(\square\)

Folgerung: Das Polynom \(c_{\acute{n}}x^{\acute{n}} + … + c_1x + c_0\) lässt sich in \(\mathcal{O}({\vartheta}^3)\) mit der Messung \((x_j, y_j) \in {}^{\omega}\mathbb{R}^2\) und den Residuen \(r_j \in {}^{\omega}\mathbb{R}\) für \(j = 1, …, m\) und \((c, x)^T \in {}^{\omega}\mathbb{R}^{\grave{n}}\) aus \(y = r + Xc\) und \(X^Tr = 0\) mit \(X = \left(x_j^k\right)\) für \(k = 0, …, \acute{n}\) bestimmen (auch unter zusätzlichen Nebenbedingungen6vgl. a. a. O., S. 196 ff.).\(\square\)

Folgerung: Die DFT (s. Wissenschaftliches Rechnen) kann die Padé-Approximation \({a(x)}/{b(x)}={\LARGE{\textbf{+}}}_{m=0}^n {a_m x^m} /\) \({\LARGE{\textbf{+}}}_{m=0}^n {b_m x^m}\) von \(f(x)\) für \(x \in {}^{\omega}\mathbb{R}\) in \(\mathcal{O}({\vartheta}^3)\) bestimmen.\(\square\)

Korollar: Intex- und zweidimensionale Bisektions- oder Newtonverfahren lösen die meisten lösbaren konvexen Programme min \(\{{f}_{1}(x) : x \in {}^{\omega}\mathbb{R}^{n}, {({f}_{2}(x), …, {f}_{m}(x))}^{T} \le 0\}\) mit konvexen Funktionen \({f}_{k} \in {}^{\omega}\mathbb{R}\) für \(k = 1, …, m\) in polynomieller Laufzeit, wenn die Anzahl der Operanden \({x}_{j}\) der \({f}_{k} \le {\omega}^{\nu-3}\) ist und ein \(x\) existiert mit \({f}_{k}(x) < 0\) für alle \(k > 1\)7vgl. Geiger, Carl; Kanzow, Christian: Theorie und Numerik restringierter Optimierungsaufgaben; 1. Aufl.; 2002; Springer; Berlin, S. 56 ff..\(\square\)

Korollar: Das LP max \(\{x : y_{\acute{m}} – xy_m = b_{\acute{m}}, y_n = y_0 = 0, x \le x_0 \in {}^{\omega}\mathbb{R}\}\) kann für \(m = 1, …, n\) nach dem Horner-Schema ein \(b_{\acute{n}}x^{\acute{n}} + … + b_1x + b_0 = 0\) lösendes \(x \in {}^{\omega}\mathbb{R}\) in \(\mathcal{O}({\vartheta}^2)\) bestimmen und abbrechende TRn erlauben jedes stetige Ungleichungssystem mit konvexer Lösungsmenge analog zu lösen.\(\square\)

© 2008-2021 by Boris Haase

Seitenbeginn

Literatur

Literatur
1 Dantzig, George B.: Lineare Programmierung und Erweiterungen; 1. Aufl.; 1966; Springer; Berlin.
2 vgl. Vanderbei, Robert J.: Linear Programming; 3rd Ed.; 2008; Springer; New York, S. 63.
3 a. a. O., S. 60 – 65.
4 Pan, Ping-Qi: Linear Programming Computation; 1st Ed.; 2014; Springer; New York, S. 580 f.
5 Vasuki, A: Nature-Inspired Optimization Algorithms; 1st Ed.; 2020; CRC Press; Boca Raton.
6 vgl. a. a. O., S. 196 ff.
7 vgl. Geiger, Carl; Kanzow, Christian: Theorie und Numerik restringierter Optimierungsaufgaben; 1. Aufl.; 2002; Springer; Berlin, S. 56 ff.