Linear Programming

Linear Programming
Linear Programming

In the following section, Set Theory, Topology and Nonstandard Analysis are presupposed.

Diameter theorem for polytopes: Every diameter of an \(n\)-dimensional polytope defined by \(m\) constraints for \(m, n \in {}^{\omega}\mathbb{N}_{\ge 2}\) is at most \(2(m + n – 3)\).

Proof: At most \(\acute{m}\) hyperplanes can be assembled into an incomplete cycle of dimension 2 and there are at most \(n – 2\) alternatives sidewards in all remaining dimensions. Overcoming every minimal distance requires at most two edges and yields the factor 2.\(\square\)

Definition: Let \(\vartheta := \omega\,{}_{\epsilon}\omega\) and \({|| \cdot ||}_{1}\) the 1-norm. If \(x \in X \subseteq {}^{\omega}\mathbb{R}^{\omega}\) passes to the next step, it is followed by \({}^{*}\) where \(\Delta x := {x}^{*} – x\). The unit vector of \(x \in X^{*}\) is \({_1}{x} := x/||x||\). Let \(x, y \in X \) in \((x, y, …)^T\) be row vectors. If a method requires computation time in seconds and memory in bits of \(\mathcal{O}({\omega}^{\mathcal{O}(1)})\), it is polynomial. If one of both quantities is \(\mathcal{O}({\epsilon}^{|\mathcal{O}(\omega)|})\), the method is exponential. Let the eigenproduct eig(\(Q\)) (yet: determinant det(\(Q\))) of a square matrix \(Q\) be the product of \(Q\)’s eigenvalues.\(\triangle\)

Theorem: The simplex method1 may solve an LP almost always in \(\mathcal{O}({\omega}^3\vartheta)\).

Proof and algorithm: Let LP max \(\{{c}^{T}x : c \in {}^{\omega}\mathbb{R}^{n}, x \in P\}\) have the feasible domain \(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}^{*}\}\). Its dual gives \({x}^{*} \in {}^{\omega}\mathbb{R}_{\ge 0}^{m}\). A (big) \(x^+ \in {}^{\omega}\mathbb{R}_{\geq0}\) achieves \(x^* \in {}^{\omega}\mathbb{R}_{\geq0}^n\) for \(x := x^* – {x^+}{\upharpoonleft}_n\). Solving max \(\{-h \in {}^{\omega}\mathbb{R}_{\le 0} : Ax – b \le h{\upharpoonleft}_m\}\) yields \(x \in P_{\ge 0}\) if \(b \ge 0\) does not hold. Let \(\left|\text{min} \left\{\complement_1^m\ b_i\right\}\right|\) be the initial value for \(h\) and 0 its target value. Start with \(x := 0\). Only one pivot step causes that \({b}^{*} \ge 0\). Let \(i, j, k \in {}^{\omega}\mathbb{N}^{*}\) and let \({a}_{i}^{T}\) the \(i\)-th row vector of \(A\). If \({c}_{j} \le 0\) for all \(j\), the LP is solved.

If for some \({c}_{j} > 0\) also \({a}_{ij} \le 0\) for all \(i\), the LP is positively unbounded. If for some \({c}_{j} = 0\) also \({a}_{ij} \le 0\) for all \(i\), drop \({c}_{j}\) and \({A}_{.j}\) as well as \({b}_{i}\) and \({a}_{i}\), but only when \({a}_{ij} < 0\) holds. The inequality \({a}_{ij}{x}_{j} \ge 0 > {b}_{i}\) for all \(j\) has no solution, too. Dividing all \({a}_{i}^{T}x \le {b}_{i}\) by \(||{a}_{i}||\) and all \({c}_{j}\) and \({a}_{ij}\) by the minimum of \(|{a}_{ij}|\) such that \({a}_{ij} \ne 0\) for each \(j\) maybe repeated and reversed later. Constraints with \({a}_{i} \le 0\) may always be removed.

For each \({c}_{j} > 0\) and non-base variable \({x}_{j}\) select \({b_k}/{a}_{kj} :=\) min \(\{{b_i}/{a}_{ij} : {a}_{ij} > 0\}\) to obtain a next \({x}^{*} \in P_{\ge 0}\) where \({x}_{j}^{*} = {x}_{j} + {b}_{k}/{a}_{kj}\). To select the steepest edge, pick the pivot \({a}_{kj}\) corresponding to \({x}_{j}\) that maximises \({c}^{T}{_1}{\Delta x}\) or \({c}_{j}^{2}/\grave{v}_{j}^{2}\) for \(v_j := ||{A}_{.j}||\) in the \(k\)-th constraint. Multiple maxima allow to use the rule of best pivot value max\({}_{k,j} {c}_{j}{b}_{k}/{a}_{kj}\) or (slower) the smallest angle min \({{_{(1)}}1{\upharpoonleft}^{T}_n}\;{_1}{c}^{*}\).

If \(c^Tx^* = c^Tx\) is repeated, relax all constraints with \({b}_{i} = 0\) by \(d \in {}^{\omega}\mathbb{R}_{>0}\). For \({b}_{i} := ||{a}_{i}||\), drop \(d\) in the tableau. If another multiple vertex is encountered for once, simply add \(||{a}_{i}||\) to every \({b}_{i}\). Leaving it may require to solve an LP with \(c > 0\) and \(b = 0\). Then the relaxation (perturbation) is reverted. This task can also finally arise, if further solutions of the LPs are to be determined for at least two \(c_j = 0\).

This avoids having to compute several minimal relaxations that differ strongly as the lexicographical method does. The rectangle rule simply allows to compute \({c}_{j}^{*}, {a}_{ij}^{*}\) and \({b}_{i}^{*}\)2. Bland’s rule is not optimal, too. No pivoting rule can prevent the simplex method from being exponential as the case may be: An exponential “drift” for e. g. Klee-Minty or Jeroslow polytopes may always provide an unfavourable edge. The result follows in accordance with the state of research.\(\square\)

Bit-square theorem (bit\({}^2\) theorem for short): Repetitions allow more than \(\mathcal{O}(n^2)\) multiplications and additions to be reduced to \(\mathcal{O}(n^2)\) for an algorithm of numbers of bit length \(m \le n\) with \(m, n \in {}^{\omega}\mathbb{N}^{*}.\square\)

Examples: The matrix product \((c_{ij}) := \left({\LARGE{\textbf{+}}}_{k=1}^n{a_{ik}b_{kj}}\right)\) of \((a_{ij}), (b_{ij}) \in {}^{\omega}\mathbb{C}^{n \times n}\) over products \({\LARGE{\textbf{+}}}_{k=1}^n{b_{jk}4^{km}}\) with \(a_{ij}\) and \(i, j \in {}^{\omega}\mathbb{N}_{\le n}^{*}\), determining an \(n\)-dimensional Taylor series (see Nonstandard Analysis), the Gauss-Jordan algorithm, the simplex and the Gram-Schmidt orthogonalisation method.

Remark: If \(m\) is too large, it will be appropriately decreased computationally. Coding quaternions enables concerted matrix multiplication in the cloud by comparing simultaneously.

Theorem: Intex methods solve every solvable LP in \(\mathcal{O}({\omega}^2\vartheta)\).

Proof and algorithm: Maybe \({b}^{T}y – {c}^{T}x \le 0, Ax \le b\) as well as \({A}^{T}y \ge c\) are normalised and scaled. By the strong duality theorem3, LP min \(\{r \in {}^{\omega}\mathbb{R}_{\ge 0} : r{\upharpoonleft}_m \ge Ax – b, r{\upharpoonleft}_n \ge c – A^Ty, r \ge b^Ty – c^Tx, b, y \in {}^{\omega}\mathbb{R}_{\ge 0}^m,\) \( c, x \in {}^{\omega}\mathbb{R}_{\ge 0}^n\} = 0\) for the residual \(r\) and the starting vector \((x^0, y^0, r^0)^T\) solves also the LPs max \(\{c^Tx : Ax \le b\}\) und min \(\{b^Ty : A^Ty \ge c\}\) in \(\mathcal{O}({}_2r^3dmn)\), if possible.

Continued interpolations using the bisection method of an inner loop approximate the respective geometric centres of the feasible regions and an extrapolation through the last two ones of the outer loop determines a new residual \(r\). If the first phase ends with \(r \ne 0\) (fixed point!), a second phase solves the LP with \(\mathcal{O}({}_2r)\) steps via two-dimensional programmes in \(\eta \in [0, 1]\) and an increased \(r\), which use \(\eta\) times \(x\) or \(y\), and remain feasible until \(r = 0.\square\)

Remark: Details are only published when no misuse for non-transparent or bad decisions must be feared. Simplex method and face algorithm4 solve LPs faster for small \(m\) and \(n\). Intex methods are faster than all known (worst-case) LP-algorithms in \(\mathcal{O}(\omega^{37/18}\vartheta)\). Optimising it for distributed computing in \({}^{\nu}\mathbb{R}^{\nu}\) reduces the runtime to \(\mathcal{O}(1)\).

Remark: Rounding errors provoked in \(b\) and \(c\) make the solution unique. Integer solutions \(v_j\) satisfying \(\acute{v}_0 \le \lfloor v_j \rfloor – v_j\) prove them successful for (mixed) integer problems and (non-) convex (Pareto) optimisation (according to nature5). Numbers of length \(\mathcal{O}({\omega})\) can only be processed in \(\mathcal{O}(\vartheta)\) as is generally known. The following results can also be easily transferred to complex numbers:

Conclusion: 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}\}\) can determine for the first solution \({x}^{o}\) any second one in \(\mathcal{O}({\vartheta}^3)\) where \({y}^{o}\) may be treated analogously.\(\square\)

Conclusion: 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*}\}\) can determine for every \(j = 1, …, n\) the eigenvalue \(\lambda_j \in {}^{\omega}\mathbb{R}\) and the eigenvector \(x_j \in {}^{\omega}\mathbb{R}^{n}\) of matrix \(A \in {}^{\omega}\mathbb{R}^{n \times n}\) in \(\mathcal{O}({\vartheta}^3).\square\)

Conclusion: LPs max \(\{{x}_{j} : Ax = 0\}\) yield all solutions to \(Ax = b\) in \(\mathcal{O}({\vartheta}^3)\). Matrix \(A\) is regular if and only if LP max \(\{{||x||}_{1} : Ax = 0\} = 0.\square\)

Conclusion: Let \({\alpha }_{j} := {A}_{.j}^{-1}\) for matrix \({A}^{-1} \in {}^{\omega}\mathbb{R}^{n \times n}\) and let \({\delta}_{ij}\) the Kronecker delta where \(j = 1, …, n\). A regular \(A\) has an eigenproduct \(\ne 0\) and allows every \({A \alpha }_{j} = {({\delta}_{1j}, …, {\delta}_{nj})}^{T}\) to be solved in \(\mathcal{O}({\vartheta}^3).\square\)

Conclusion: The polynomial \(c_{\acute{n}}x^{\acute{n}} + … + c_1x + c_0\) can be determined in \(\mathcal{O}({\vartheta}^3)\) by measurement \((x_j, y_j) \in {}^{\omega}\mathbb{R}^2\) and residues \(r_j \in {}^{\omega}\mathbb{R}\) for \(j = 1, …, m\) and \((c, x)^T \in {}^{\omega}\mathbb{R}^{\grave{n}}\) from \(y = r + Xc\) and \(X^Tr = 0\) where \(X = \left(x_j^k\right)\) for \(k = 0, …, \acute{n}\) (also under additional constraints6).\(\square\)

Conclusion: The discrete Fourier transform (see Scientific Computing) can determine the Padé approximant \({a(x)}/{b(x)}={{\LARGE{\textbf{+}}}_{m=0}^n {a_m x^m}}/{{\LARGE{\textbf{+}}}_{m=0}^n {b_m x^m}}\) of \(f(x)\) for \(x \in {}^{\omega}\mathbb{R}\) in \(\mathcal{O}({\vartheta}^3).\square\)

Corollary: Intex methods and two-dimensional bisection or Newton’s methods may solve most solvable convex programmes min \(\{{f}_{1}(x) : x \in {}^{\omega}\mathbb{R}^{n}, {({f}_{2}(x), …, {f}_{m}(x))}^{T} \le 0\}\) where all \({f}_{k} \in {}^{\omega}\mathbb{R}\) are convex functions for \(k = 1, …, m\) in polynomial runtime, if the number of operands \({x}_{j}\) of all \({f}_{k}\) is \(\le {\omega}^{\nu-3}\) and if an \(x\) exists such that \({f}_{k}(x) < 0\) for all \(k > 1\)7.\(\square\)

Corollary: LP max \(\{x : y_{\acute{m}} – xy_m = b_{\acute{m}}, y_n = y_0 = 0, x \le x_0 \in {}^{\omega}\mathbb{R}\}\) can determine for \(m = 1, …, n\) by Horner scheme an \(x \in {}^{\omega}\mathbb{R}\) solving \(n\)-polynomial \(b_{\acute{n}}x^{\acute{n}} + … + b_1x + b_0 = 0\) in \(\mathcal{O}({\vartheta}^2)\). Aborting TS can analogously solve every continuous inequation system with convex solution set.\(\square\)

Theorem for the Strassen algorithm: Computing the matrix product \(AA^T\) shortens the original runtime \(T(n) = \mathcal{O}(n^{(_2 7)})\) roughly by \(\tilde{3}\) for sufficiently big \(n := 2^k, k \in \mathbb{N}^*\) and the matrix \(A \in \mathbb{C}^{n \times n}\) due to the GS and\[A := \begin{pmatrix} A_{11} & A_{12} \\ A_{21} & A_{22} \end{pmatrix} \text{ as well as } AA^T = \begin{pmatrix} A_{11}A_{11}^T+A_{12}A_{12}^T & A_{11}A_{21}^T+A_{12}A_{22}^T \\ A_{21}A_{11}^T+A_{22}A_{12}^T & A_{21}A_{21}^T+A_{22}A_{22}^T \end{pmatrix}.\square\]

code of the simplex method

© 2008-2024 by Boris Haase

top

References

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