1296.3 Page 1296.3:Algorithms in the Real WorldLinear and Integer Programming II– Ellipsoid algorithm– Interior point methods296.3 Page 2Ellipsoid AlgorithmFirst known-to-be polynomial-time algorithm for linear programming (Khachian 79)Solvesfind xsubject to Ax ≤ bi.e., find a feasible solutionRun Time:O(n4L), where L = #bits to represent A and bProblem in practice: always takes this much time.296.3 Page 3Reduction from general caseTo solve:maximize: z = cTxsubject to: Ax ≤ b, x ≥ 0Convert to:find: x, y subject to: Ax ≤ b,-x ≤ 0-ATy ≤ –c-y ≤ 0-cTx +yTb ≤ 0Could add constraint-cTx ≤ -z0, do binarysearch over various values of z0to find feasible solution with maximum z, but approach based on dual gives direct solution. 296.3 Page 4Ellipsoid AlgorithmConsider a sequence of smaller and smaller ellipsoids each with the feasible region inside.For iteration k:ck= center of EkEventually ckhas to be inside of F, and we are done.ckFFeasible region2296.3 Page 5Ellipsoid Algorithmfind smallest ellipsoid that includes constraintTo find the next smaller ellipsoid:find most violated constraint akckFFeasible regionak)12/(1121)()(++=nkkEVolEVol296.3 Page 6Interior Point MethodsTravel through the interior with a combination of1. An optimization term(moves toward objective)2. A centering term(keeps away from boundary)Used since 50s for nonlinear programming.Karmakar proved a variant is polynomial time in 1984x1x2296.3 Page 7MethodsAffine scaling: simplest, but no known time boundsPotential reduction: O(nL) iterations Central trajectory: O(n1/2L) iterationsThe time for each iteration involves solving a linear system so it takes polynomial time. The “real world” time depends heavily on the matrix structure.296.3 Page 8Example timesCentral trajectory method (Lustic, Marsten, Shanno 94)Time depends on Cholesky non-zeros (i.e., the “fill”)6.7M.2M.3M1.2MCholeskynon-zeros92526457712364time (sec)58536466iterations80K183K189K186Knon-zeros19x12K43x107K9x57K13x31Ksize (K)initialcarcontinentfuel3296.3 Page 9AssumptionsWe are trying to solve the problem:minimize z = cTxsubject to Ax = bx ≥ 0296.3 Page 10Outline1. Centering Methods Overview2. Picking a direction to move toward the optimal3. Staying on the Ax = b hyperplane (projection)4. General method5. Example: Affine scaling6. Example: potential reduction7. Example: log barrier296.3 Page 11Centering: option 1The “analytical center”:Minimize: y = -Σi=1nlg xiy goes to infinity as x approaches any boundary.x1x2x4x3x5296.3 Page 12Centering: option 2Elliptical Scaling:(c1,c2)122222121=+cxcxDikin EllipsoidThe idea is to bias spaced based on the ellipsoid.More on this later.4296.3 Page 13Finding the Optimal solutionLet’s say f(x) is the combination of the “centering term” c(x) and the “optimization term” z(x) = cTx. We would like this to have the same minimum over the feasible region as z(x) but can otherwise be quite different.In particular c(x) and hence f(x) need not be linear.Goal: find the minimum of f(x) over the feasible region starting at some interior point x0Can do this by taking a sequence of steps toward the minimum.How do we pick a direction for a step?296.3 Page 14Picking a direction: steepest descentOption 1: Find the steepest descent on x at x0by taking the gradient:Problem: the gradient might be changing rapidly, so local steepest descent might not give us a good direction.Any ideas for better selection of a direction? )(0xf∇296.3 Page 15Picking a direction: Newton’s methodTo find the minimum of f(x) take the derivative and set to 0.200000))((''21))((')()( xxxfxxxfxfxf −+−+≈Consider the truncated Taylor series:))(('')('0000xxxfxf −+=)('')('000xfxfxx −=In matrix form, for arbitrary dimension:TxfxFxx )())((10∇−=−)()( xfxF ∇×∇=Hessian 296.3 Page 16Next Step?Now that we have a direction, what do we do?5296.3 Page 17Remaining on the support planeConstraint: Ax = b A is a n x (n + m) matrix.The equation describes an m dimensional hyperplanein a n+m dimensional space.The hyperplane basis is the null space of AA = defines the “slope”b = defines an “offset”x2x13x1+ 2x2= 3x1+ 2x2= 44296.3 Page 18ProjectionNeed to project our direction onto the plane defined by the null space of A.PccAAAAIAcAAAcncdTTTT=−=−=−=−−))(()(11()matrix" projection" the1=−=−AAAAIPTTWe want to calculate Pccnd296.3 Page 19Calculating PcPc = (I – AT(AAT)-1A)c = c – ATwwhere ATw = AT(AAT)-1Ac giving AATw = AAT(AAT)-1Ac = Acso all we need to do is solve for w in: AATw = AcThis can be solved with a sparse solver as described in the graph separator lectures.This is the workhorse of the interior-point methods.Note that AAT will be more dense than A.296.3 Page 20Next step?We now have a direction c and its projection d onto the constraint plane defined by Ax = b.What do we do now?To decide how far to go we can find the minimum of f(x) along the line defined by d. Not too hard if f(x) is reasonably nice (e.g., has one minimum along the line).Alternatively we can go some fraction of the way to the boundary (e.g., 90%)6296.3 Page 21General Interior Point MethodPick start x0Factor AAT(i.e., find LU decomposition)Repeat until done (within some threshold)– decide on function to optimize f(x)(might be the same for all iterations)– select direction d based on f(x)(e.g., with Newton’s method)– project d onto null space of A (using factored AAT and solving a linear system)– decide how far to go along that directionCaveat: every method is slightly different296.3 Page 22Affine Scaling MethodA biased steepest descent.On each iteration solve:minimize cTysubject to Ay = 0yTD-2y ≤ 1Note that: 1. Ellipsoid is centered around current solution x2. y is in the null space of A and can therefore be used as the direction d without projection3. we are optimizing in the desired direction cTWhat does the Dikin Ellipsoid do for us?=O321000000xxxDDikin ellipsoid296.3 Page 23Dikin EllipsoidFor x > 0 (a true “interior point”),=O321000000xxxD122222221212≤+++=−nnTxyxyxyyDy Lx1x2x4x3x5This constraint prevents y from crossing any face. Ay=0 keeps y on the right hyperplane.Optimal value on boundary of ellipsoid due to convexity.Ellipsoid biases searchaway from corners.yx296.3 Page 24Finding Optimalminimize cTysubject to Ay = 0yTD-2y ≤ 1But this looks harder to solve than a linear