15.093J Optimization Methods Lecture 22: Barrier Interior Point Algorithms1 Outline Slide 1 1. Barrier Methods 2. The Central Path 3. Approximating the Central Path 4. The Primal Barrier Algorithm 5. The Primal-Dual Bar rier Algorithm 6. Computational Aspects of IPMs 2 Barrier methods Slide 2 min f(x) s.t. gj(x) ≤ 0, j = 1, . . . , p hi(x) = 0, i = 1, . . . , m S = {x| gj(x) < 0, j = 1, . . . , p, hi(x) = 0, i = 1, . . . , m} 2.1 Strategy Slide 3 • A barrier function G(x) is a continous function with the property that is approaches ∞ as one of gj(x) approaches 0 from negative values. • Examples: p p �� 1 G(x) = − log(−gj(x)), G(x) = − gj(x)j=1 j=1 Slide 4 Consider a sequence of µk: 0 < µk+1 < µk and µk 0. • →• Consider the problem x k = argminx∈S � f(x) + µ kG(x) � • Theorem If Every limit point xk generated by a barrier method is a global minimum of the original constrained problem. 1� � x* . . . . x(0.01) x(0.1) x(1) x(10) . central path . analytic center c 2.2 Primal path-following IPMs for LO Slide 5 (P ) min c�x (D) max b�p s.t. Ax = b s.t. A�p + s = c x ≥ 0 s ≥ 0 Barrier problem: n min Bµ(x) = c � x − µ log xj j=1 s.t. Ax = b Minimizer: x(µ) 3 Central Path Slide 6 • As µ varies, minimizer s x(µ) form the c e ntral path ∗ • limµ→0 x(µ) exists and is an o ptimal solution x to the initial LP • For µ = ∞, x(∞) is called the analytic center n min log xj− j=1 s.t. Ax = b Slide 7 2� � � .. x1 x2 x3 Q the central path the analytic center of Q the analytic center of P (1/3, 1/3, 1/3) (1/2, 0, 1/2) P • Q = x | x = (x1, 0, x3), x1 + x3 = 1, x ≥ 0}, set of optimal solutions to original LP • The analytic center of Q is (1/2, 0, 1/2) 3.1 Example Slide 8 min x2 s.t. x1 + x2 + x3 = 1 x1, x2, x3 ≥ 0 min x2 − µ log x1 − µ log x2 − µ log x3 s.t. x1 + x2 + x3 = 1 min x2 − µ log x1 − µ log x2 − µ log(1 − x1 − x2). x1(µ) =1 −x2(µ) 2 1 + 3µ − 1 + 9µ2 + 2µ x2(µ) = 2 x3(µ) =1 −x2(µ) 2 The analytic center: (1/3, 1/3, 1/3) Slide 9 3.2 Solution of Central Path Slide 10 • Barrier problem for dual: n max p �b + µ log sj j=1 s.t. p �A + s � = c � 3� � • Solution (KKT): Ax(µ) = b x(µ) ≥ 0 A� p(µ) + s(µ) = c s(µ) ≥ 0 X(µ)S(µ)e = eµ Slide 11 Theorem: If x ∗ , p ∗, and s ∗ satisfy optimality conditions, then they are • optimal solutions to problems primal and dual barrier problems. • Goal: Solve barrier problem n min Bµ(x) = c�x − µ log xj j=1 s.t. Ax = b 4 Approximating the central path Slide 12 ∂Bµ(x) µ ∂xi = ci −xi ∂2Bµ(x) µ = x∂x2 i i 2 ∂2Bµ(x) = 0, i = j∂xi∂xj �Given a vector x > 0: Slide 13 n � ∂Bµ(x)Bµ(x + d) ≈ Bµ(x) + di + ∂xii=1 1 n ∂2Bµ(x)didj2 ∂xi∂xji,j=1 = Bµ(x) + (c � − µe �X−1)d + 12µd�X−2d X = diag(x1, . . . , xn) Slide 14 Approximating problem: min (c � − µe �X−1)d +1 µd�X−2d 2s.t. Ad = 0 Solution (from Lagrange): c −µX−1 e + µX−2d −A� p = 0 Ad = 0 Slide 15 4� � � � � � � � x* . . . . x(0.01) x(0.1) x(1) x(10) . central path . analytic center c • System of m + n linear equations, with m + n unknowns (dj, j = 1, . . . , n, and pi, i = 1 , . . . , m). Solution: • d(µ) = I − X2A�(AX2A�)−1A xe −µ 1 X2 c p(µ) = (AX2A�)−1A(X2 c − µxe) 4.1 The Newton connection Slide 16 • d(µ) is the Newton direction; process of calculating this direction is called a Newton step • Starting with x, the new primal solution is x + d(µ) The corresponding dual solution becomes (p, s) = � p(µ), c −A� p(µ) � • • We then decrease µ to µ = αµ, 0 < α < 1 4.2 Geometric Interpretation Slide 17 • Take one Newton step so that x would be close to x(µ) Measure of closeness �� �� • �� 1 �� �� µ XSe − e�� ≤ β, 0 < β < 1, X = diag(x1, . . . , xn) S = diag(s1, . . . , sn) As µ 0, the complementarity slackness condition will be satisfied • →Slide 18 5� � � � � � 5 The Primal Barrier Algorithm Slide 19 Input (a) (A, b, c); A has full row rank; (b) x0 > 0, s0 > 0, p0; (c) optimality tolerance � > 0; (d) µ0, and α, where 0 < α < 1. Slide 20 1. (Initialization) Start with some primal and dual feasible x0 > 0, s0 > 0, p0, and set k = 0. 2. (Optimality test) If (sk)�xk < � stop; else go to Step 3. 3. Let Xk = diag(x1 k , . . . , x nk), µ k+1 = αµk Slide 21 4. (Computation of directions) Solve the linear s ystem µk+1X−k 2d − A� p = µk+1X−k 1 e −c Ad = 0 5. (Update of solutions) Let x k+1 = x k + d, k+1 p = p, s k+1 = c −A� p. 6. Let k := k + 1 and go to Step 2. 5.1 Correctness Theorem Given α = 1 −√√ββ + −√βn, β < 1, (x0 , s0 , p0), (x0 > 0, s0 > 0): Slide 22 �� 1 �� �� X0S0e − e�� ≤ β. �� µ0 �� Then, after √β + √n (s0)�x0(1 + β)K = log √β −β�(1 −β) iterations, (xK , sK , pK) is found: (s K)� x K ≤ �. 6�� � � � � � � � � 5.2 Complexity Slide 23 • Work per iteration involves solving a linear system with m + n equations in m + n unknowns. Given that m ≤ n, the work per iteration is O(n3). �0 = (s0)�x0: initial duality gap. Algorithm needs • O √n log �0 iterations to reduce the duality gap from �0 to �, with O(n3) arithmetic operations per iteration. 6 The Primal-Dual Barrier Algorithm 6.1 Optimality Conditions Slide 24 Ax(µ) = b x(µ) ≥ 0 A�p(µ) + s(µ) = c s(µ) ≥ 0 sj(µ)xj(µ) = µ or X(µ)S(µ)e = eµ X(µ) = diag x1(µ), …