Master`s Thesis Implementation of a bundle algorithm
Contents
1. 6 2 3 3 The primal solution 7 Dual Algorithms 8 3 1 Subgradient methods o 9 3 2 Conjugated subgradient methods 9 3 3 e subgradient methods o o 10 3 4 Cutting plane methods o e 11 3 5 Stabilized cutting plane methods 12 3 6 Method connection scheme 2 0004 13 Bundle Algorithms 14 Aol TESTA GEBIES 22 0 ca O a Go A lt E eis 15 4 2 Performance results eee ee eee 16 Manual 17 5il Preparations zy ies Ge ow ww ee Re od a te oe kh 17 Did TRUE coves eae kde He Soh A AANA ae ee eee lee Y 17 Deeds OTAC eas sha atk Sa ek Bo Re fee eae ee a ace i 18 9 2 2 Parameters 2 24444 228 aa ee beeen Greed 19 5 2 3 Main anda IEA AR BE eS 19 52A Makefile car E E tele et Be a 19 5 3 EXamples Lia Ol as eee ee dedo ed 20 Drow Example ana ma tt A A eth ey ad 20 53 2 Example 4 4 2 2 a et lay de tada 23 ody Tutorials 9 4 i amp c4 A RIA AA ES at ak 25 5 4 1 Set Covering Problem SCP oooo o o oo 25 5 4 2 Generalized Assignment Problem GAP 26 5 5 Subgradient implementation sosoo e 27 A The dual problem when there are primal inequality constraints 28 Optimality conditions when there are primal inequality con straints 29 Implementation files and their connections 30 More Set Covering Problem SCP results 32 References 33 1 Introduction In this2. Definition 1 Filling Property The filling property for 2 4 is said to hold at u E R if 00 u is the convex hull of the set G u defined in 7 So when does the filling property hold The question is rather technical and will not be answered in this report But we will list some situations when it does hold To get a positive answer some extra requirements of a topological nature must be added to the primal problem 2 Theorem 3 The filling property of Definition 1 holds at any u R when X is a compact set on which f and each c are continuous in particular when X is a finite set such as in the usual Lagrangian relaxation of combinatorial problems in linear programming and in quadratic programming more generally in problems where f and each cj are quadratic or even lp norms 1 lt p lt 00 When the filling property holds at u any subgradient of 9 at u can be written as a convex combination of constraint values at sufficiently but finitely many x s in X u For example g E Ou g 5 arg c cz k where the ay s are convex mulitpliers ic X pap 1 and a gt 0 Vk A u minimizing 0 5 has the property 0 00 u which means that there exists sufficiently many x s in X u such that the c x s have 0 in their convex hull Proposition 1 Let u solve the dual problem 5 and suppose the filling prop erty holds at u Then there are at most m 1 points zg and convex
3. Definition 2 Bundle The bundle B is a container of information about the function A ui g ui Eli TE where I is the number of pairs 0 gi in the bundle Notice that different methods store different numbers of pairs 0i gi thus I can be the number of iterations made or a lower number i e we don t save all pairs 0 g9 that have been generated There are essentially two types of methods subgradient and cutting plane methods lying at opposite sides of a continuum of methods In between we find among others bundle methods We will first describe methods derived from subgradients and then from cutting planes Last we will present a scheme visualizing the connections between the methods In the following we will denote by gx the vector g uz obtained from the oracle at some iterate uz Further we let ux be the current iterate 3 1 Subgradient methods Subgradient methods are very common because of their simplicity At iteration K select tg gt 0 and set UK 1 UK tKgK No control is made on gx to check if it is a good direction or even a descent direction The only control of the method s behaviour is to choose the stepsize tx In our implementation of a subgradient method in Chapter 5 5 we use the following step size rule M AA b lt M k Lenght b k btk gt 0 O lt uws lt 00 Eddy pi tk which is proved to converge in the sense that 0 ux 0 when K oo The pr
4. The base for this technique is that the steepest descent direction with respect to 0 at u is the minimum norm vector g 00 u the problem 9 is therefore an approximation of the problem of finding the steepest descent direction A major drawback with this method is that every time a good improvement is possible and the current point u is moved 4 must be emptied loosing all the information about 0 This is because the elements of 3 do not really need to be subgradients of the new current point More information can be found in Frangioni 3 pp 5 6 3 3 e subgradient methods Instead of using the somewhat nontheoretical rules by which subgradients are used in the conjugated subgradient methods we can instead define approximate subgradients as follows Definition 3 e subgradient Lete gt 0 The vector g is an e subgradient of 0 atu if 0v gt O u g v u e for any ve R 90 Figure 2 e subgradient g and its linearization error a Now we fix some e gt 0 and solve 9 further restricted to the known subgradients at the current point u min 191 g E conv gi ai lt i E BP where a is the linearization error see Figure 2 Instead of eliminating all the gis with a gt e the a can be used as weights of the corresponding vector In that way any g with a big a will have a small influence on the direction min gyll X api lt e p 9 10 1 B 1
5. w whose effect is to pull ux 41 toward u thereby stabilizing the algorithm We then let t gt 0 dynamically if we want decide the effect of the stabilizing term 1 O u u zl all for all u R 14 Minimizing 6 over R is equivalent to the quadratic problem 1 min r zl all u r RH 15 r gt lui 9 u u i for allie 8 which always has a unique optimal solution ux 1 7k 1 The two positive numbers 0 t O ux41 G y 1 8 blur 8 lla l still represent the expected decrease 0 U 0 uK 1 To complete the iteration we have to update the stability center u This will only be done if the improvement is large enough relative to 6 O ux41 lt O a Kd r EJO 11 16 k being a fixed tolerance The bigger K the better must the model be to satisfy the condition if satisfied Ux 1 If the condition is not satisfied we will at least work with a better model in the next iteration This is the key concept in bundle methods The bundle is updated after each iteration The new element Ox 1 9K 1 from the oracle is added while old information is deleted if not in use for a certain time For example if an element 0 9 8 from a previous iteration does not define an active piece of the model defined by 6 during many iterations in other words strict inequality holds in the corresponding constraint in 15 for several consecutive
6. 6 and 7 The oracle is a C class made up by two files Oracle h and Oracle C Oracle h contains the class declarations and Oracle C the class implementa tions Here follows a compressed Oracle h 01 class Oracle 02 03 private Private 04 IloEnv env Global Concert Technology variables 05 IloModel model 06 07 TloInt uSize Global necessary variables 08 IloInt xDims 09 E 10 public Public TIT 11 Oracle char filename 12 Oracle 13 void Solve const double uVector double amp fValue 14 double sgVector double xVector 15 16 The problem is modelled in the constructor Oracle row 11 and solved in the function Solve row 13 The communication between them is handled by the global variables in the private part of the declaration Because the constructor is only run once and Solve is run every time the algorithm asks for a new piece to the bundle we try to do as much work as possible in the constructor such as reading data and creating the model Observe that the constructor takes one argument a filename which can be used if you build a solver for a class of problems with data read from different files The function Solve takes a u as in parameter and returns the function value 0 u the subgradient g u and the primal x point of the solution to the inner problem in 4 The x point is necessary to build a primal solution when the dual is found as in 8 18 The deta
7. If we compare these results with them from the bundle implementation we see that the optimal solution is never reached even after 1000 iterations The subgradient implementation has also large problems to adapt to different scales of problems notice for example the SCP small txt With these result as proof we can truly say that the bundle implementation is far superior to this simple subgradient algorithm 27 A The dual problem when there are primal in equality constraints Theorem 4 Suppose the primal problem 2 has inequality constraints max f z st 2 X CR cls lt 0 gj 1 m Then the dual variable u must be positive and the dual problem becomes min 0 u 18 wER where 0 u matecx f x u c x Proof To get the form 2 introduce slack variable s and use the flexibility of Lagrangian relaxation set X X x RY so that we have to dualize max f x st TEXCR s gt 0E R c z s 0 ER The Lagrangian is L a s u f x u c a s L z u u s where L x u is the ordinary Lagrangian as though we had equalities The resulting dual function is clearly ae oo if some uj lt 0 a 0 u otherwise where is the ordinary dual function In a word the dual problem becomes 18 28 B Optimality conditions when there are primal inequality constraints Theorem 5 Suppose the primal problem 2 has inequality constraints as in Theorem
8. PB This is the first method referred to as a bundle method The size of the bundle is much larger compared to the subgradient methods above The critical part of this approach is the introduction of the parameter e which has to be dynamically updated to achieve good performance This has turned out to be difficult and it is not clear what sort of rules that should be used For example this further smoothing has been introduced oa min tl X gill opi p Ed 11 1 8B 1 8 10 which can be viewed as a Lagrangian relaxation of 10 with some addi tional smoothing where 1 t is a parameter used to dualize the constraint AE papi SE Now we have t instead as a parameter which has turned out to be a little eas ier to update than e Different rules for updating t are explained in Chapter 4 1 More information can be found in Frangioni 3 pp 6 8 3 4 Cutting plane methods In this section we will present another type of methods in which a piecewise affine model of 0 will be minimized at each iteration Every result from the oracle with uz as input defines an affine function approximating 0 from below O u gt O ux g u uz forall u R where g 00 u with equality for u uz After K iterations we can construct the following piecewise affine function 9 which underestimates 6 9 u gt Ou max 6 u g u us Bj forall u R The problem is now to calculate an optimal solutio
9. ua 4 21 213 xa 323 lt 6 1 2 3 E 1 10 ui ER we RY Dual function 4 O u min L x u u ER uzER x st 2 3a3 lt 6 21 122 123 1 10 The code Only some differences from Example 1 will be highlighted the rest follow the same structure as before When we define the primal x variables we can set the bounds directly instead of giving them as constraints to the model x primal variables x add IloIntVarArray env n 0 10 We have also defined them as integers by using the IloIntVarArray In this example we have both equality and inequality relaxed constraints i e uy is free unconstrained and uz must be positive constrained The code for the ucVector becomes Add all unconstrained multipliers u and with an InINF ucVector new unsigned int 2 ucVector 0 0 ucVector 1 InINF gt We add one element to the ucVector the first one with the name O and finish off with the end mark InINF 23 The bundle algorithm always expects a max problem to be solved but in the current example we have a min problem The solution to this is to change the sign of the function value and the subgradient when returned in the Solve function fValue with because min problem fValue cplex getValue obj sgVector with because min problem sgVector 0 cplex getValue lag1 sgVector 1 cplex getValue lag2 Results Optimal soluti
10. 4 Then the property 0 00 u is not a necessary condition for opti mality It is instead enough if the complementary conditions u gu 0 Ju ER weR where gu 00 u 19 are fulfilled which means that either uj or gu has to be 0 for all j 1 m Proof We will show the condition for one dimension then it s easy to convince us of the general case Our problem is then to show the optimality conditions when minimizing the convex 1 dimensional function min 0 u s t u gt 0 We get two scenarios Figure 5 Two scenarios for optimality From the figures we can see that we have a minimum when either gu 0 or u 0 in combination with a gu gt 0 The optimality conditions for one dimension becomes u gu 0 gu 0 These conditions must be true for every dimension and the general case becomes 19 29 C Implementation files and their connections There are three types of C files in the implementation h header files C implementation files and o compiled object files In some classes is only the object file included and not the implementation file the reason for this are restrictions in code distribution Here follows a list of all the files in the bundle implementation OPTtypes h Help class with standard types and constants OPTvect h Help class with scalar and vector functions CMinQuad h o Implementation of the TT Tall and Thin algorithm for solving the typical Qua
11. enough improvement has been obtained then 06 u u SS else it is a NS 07 update 08 update t 09 until some stopping condition Explanations to some of the rows 01 is usually set to all zeros and a call to the oracle gives the first piece to PP Other parameters control strategies for and t as well as decide the precision in the stopping test 03 This is the same as solving one of the two quadratic problems 11 with a line search or 15 A fast QP solver is therefore a must in a good bundle method implementation 04 A call to the user created oracle gives new information about at u 05 This can be a test like 16 07 8 is updated with the new piece of information from the oracle Old information is deleted if not in use for a certain time 08 t is updated increased or decreased depending on the model s accuracy Further details are found in the next subsection t strategies 09 This can be a test like 17 A more detailed version can be found in Frangioni 3 pp 63 70 14 4 1 t strategies What does the parameter t affect By studying the subproblem of the stabilized cutting plane method 15 we see that controls the strength of the quadratic term whose effect is to pull uxy1 toward u A large t will result in a weak stabilization and will be used when the model is good and a small will result in a strong stabilization and will be used when the model is bad This is w
12. report we want to solve problems of the type min 0 u we 1 where 0 R gt R is a convez and finite everywhere nondifferentiable function and U isa nonempty convex subset of R One of the most important classes of these problems is Lagrangian duals there is given implicitly by a maximization problem where u is a parameter called Lagrangian multiplier We will assume in this report that we have an oracle which returns 0 u and one subgradient glu of O at u for any given u U see Figure 1 Because of the strong connection with Lagrangian duals algorithms solving 1 are often referred to as dual algorithms In bundle methods one stores information 6 u g u about previous it erations in a set 8 the bundle to be able to choose a good descent direction d along which the next points are generated Having access to this memory makes it possible to construct a better search direction than so called subgradient methods that only utilize the information from the present point In Chapter 2 we will describe Lagrangian relaxation and see how it gives rise to problems of the type 1 In Chapter 3 we describe different types of dual algorithms including bundle methods that solve 1 Chapter 4 will be devoted to a more detailed study of the bundle algorithm Chapter 5 will contain the manual for our bundle implementation 90 real function Figure 1 The function u is defined as the maximization of a set
13. Masters Thesis Implementation of a bundle algorithm for convex optimization Reine Saljo May 11 2004 _Mathematical Sciences GOTEBORG UNIVERSITY Abstract The aim of this thesis is to implement a bundle algorithm for convex optimiza tion We will describe the technique the theory behind it and why it has been developed Its main application the solution of Lagrangian duals will be the foundation from which we will study the technique Lagrangian relaxation will therefore also be described We will further in a manual describe the implementation of the code and sim ple examples will show how to use it Tutorials have been developed for solving two problem types Set Covering Problem SCP and Generalized Assignment Problem GAP A simple subgradient algorithm has also been developed to show the benefits of the more advanced bundle technique Acknowledgements I would like to thank Professor Antonio Frangioni at Universita di Pisa Genova Udine for sharing his bundle implementation developed during his research I will also thank Professor Michael Patriksson my tutor and examiner during this Master s Thesis Contents 1 2 Introduction 3 Lagrangian Relaxation 4 2 1 General technique o 4 2 2 Minimizing the dual function 5 2 3 Primal dual relations o 00004 6 2 3 1 Everett s Theorem 6 2 3 2 Differential study of the dual function
14. ality the same Subgradient method y Conjugated subgradient method e subgradient method Lagrangian relaxation dual Bundle methods Stabilized cutting plane method 1 Cutting plane method Figure 3 Method connection scheme where 9 10 etc represent the subproblems in the stated methods in the figure 13 4 Bundle Algorithms In this section we will take a closer look at the bundle algorithm based on 11 or 15 the two dual subproblems in the bundle methods A subsection will also be devoted to t the single most important parameter for the behaviour and performance of the bundle algorithm As described earlier bundle algorithms collect information 9 u g u about previous iterations in a set P the bundle The information is given by a user created oracle the only information available about the problem This information is then used to compute a tentative descent direction d along which the next points are generated Since d need not to be a direction of descent bundle methods use some kind of control system to decide to either move the current point u SS a Serious Step or use the new information to enlarge 3 to be able to compute a better d at the next iteration NS a Null Step Here follows a generic bundle algorithm 01 initiate a 8 t and other parameters 02 repeat 03 find new direction d generate a new trial point u u 7d 04 oracle u gt 0 u g u 05 if a large
15. ax L x u ueR st 0 lt x lt 4 The code First we have to create the Concert Technology environment an instance of IloEnv to handle memory management for all Concert Technology objects IloEnv env This variable has to be available in all functions and is therefore declared in the private part of the Oracle h file Now we can start with the model located in the constructor First the basic variables m rows n cols m 1 n 2 x primal variables x add IloNumVarArray env n 0 IloInfinity 20 u dual multipliers u add IloNumVarArray env m All are global variables declared in the Oracle h file The elements m and n are the size of the x matrix used when creating arrays and in for loops The elements x and u are declared as empty arrays and to initiate them we add n elements to x and m elements to u Now let s create the model model primal object cx x 0 2 x 1 model lagrange relaxed constraint subgradient lag x 0 4 x 1 8 model maximize obj IloMaximize env cx u 0 lag model add obj model x constraints for j 0 j lt n j model add x j gt 0 model add x j lt 4 cplex extract model cplex extract model Also here is all the variables global declared in the Oracle h file The ex pressions cx and lag are instances of IloExpr holding mathematical expres sions The objective an instance of I1o0bject
16. cess ilm Observe that it contains the path to the CPLEX installation which may have changed This row can be copied and pasted from the top of the Makefile file 5 2 Files There are many files in the implementation but there are mainly four files the user will have to interact with these are Oracle C The only file containing information about the problem Oracle h Oracle header file Parameters Contains initial parameters Main C Controls among other things the printouts For a description of all the files and their connections see Appendix C 17 5 2 1 Oracle As described earlier the oracle returns a function value 0 u and a subgradient g u for any u R We have also seen that when solving the easy dual function 4 we get this information for free But still the problem has to be solved and for this task we have chosen ILOG CPLEX This gives us the obvious limitation the maximization problem defining the dual function 4 must be solvable in CPLEX i e it must be a linear problem LP or a quadratic problem QP or a linear integer problem LIP More information is found in the CPLEX manual 5 The technique we use to represent problems is called ILOG Concert Tech nology and is a C library of classes and functions making it easy to input various problems We will describe this technique through some guided exam ples in Section 5 3 but for more details you have to read the CPLEX manuals
17. d the dual problem becomes min u weER To see why read Appendix A R changes to R in all the theory above We can now draw the conclusion that the dual problem 5 fulfil the require ments of problem 1 2 2 Minimizing the dual function In this section we will state a very important theorem when minimizing the dual function The fundamental result is Theorem 1 Whatever the data f X c in 2 can be the function 0 is always conver and lower semicontinuous 1 Besides if for given u R 4 has an optimal solution x not necessarily unique then gu C is a subgradient of 0 at u bv gt O u gi v u for any v R 6 which is written as gy 00 u A consequence of this theorem is that the dual problem 5 is well posed minimizing a convex function We also see that every time we maximize the Lagrangian 4 we get for free the value 0 u and the value c 2 of a subgradient both important to solve the dual problem Observe that c x is the vector of partial derivates of L with respect to u ie VuL t u As a result solving the dual problem 5 is exactly equivalent to minimizing a convex function with the help of an oracle the dual function 4 providing function and subgradient values TA function 0 is lower semicontinuous if for any u R the smallest value of 0 v when v gt wis at least 0 u liminf f y gt f a you Lower sem
18. dratic Problem BMinQuad h o Implementation of the BTT algorithm for solving box constrained Quadratic Problems Bundle h o The bundle algorithm QPBundle h o Implementation of the QP direction finding subprob lem for the Bundle class using the specialized QP solvers implemented in the CMinQuad class and BMinQuad class Solver h C Inherits the QPBundle class and connects to the Oracle class Oracle h C The user provided problem we want to solve Uses ILOG CPLEX Main C The main program creating instances of the Solver class and the Oracle class Parameters Contains initial parameters MakeFile Contains all the compiler directives and links all the files to gether with the ILOG CPLEX environment 30 Used by all OPTtypes h OPTvect h Bundle solver structure CMinQuad BMinQuad QPBundle Main program structure Figure 6 Connection scheme for bundle implementation 31 D More Set Covering Problem SCP results Here follow some more results produced with the SCP tutorial Five prob lems with increasing size have been solved to optimality Each problem has been solved with three different strategies pure heuristic heuristic with soft long term strategy and heuristic with hard long term strategy all described in Section 4 1 All the example files can be found in lt Examples directory gt SCP More data files can be found at 8 The parameter controlling the strategy of the solver is f
19. e to do is to update the model Update u u clear for i 0 i lt m i u add uVector i Update model model remove obj obj IloMaximize env cx u 0 lag model add obj u is updated and the old objective is removed The same code as above creates the new objective and is added to the model The cplex object is automatically updated with the new model and we can now solve it Solve model cplex solve The last thing to do is to create the return values of the function Solve the objective value 6 u the subgradient g u and the x point fValue fValue cplex getValue obj sgVector sgVector 0 cplex getValue lag xVector for j 0 j lt n j xVector j cplex getValue x j Here we use the function getValue not only to read the objective but also to read the current values of different parts of the model like lag and the x point The rest of the code is left to the user to study It contains some code for printouts and some necessary get functions for the class Results Optimal solution in 5 iterations f 6 4 1 u 0 5 22 5 3 2 Example 2 The code to this example can be found in lt Examples directory gt EX2 The problem Primal problem 2 min 3x1 5x9 423 st 27 23 6 relax 1 2 2 gt 4 relax t 323 lt 6 z1 02 23 1 10 Lagrangian 3 L x u 311 5 2 4x3 u 211 3 6
20. elative increments are measured and the absolute magnitude of the obtainable increments is disregarded After some initial iterations t starts to rapidly decrease and soon the algorithm is working with a very small t resulting in extremely slow convergence The conditions for t to increase are almost never satisfied Long term strategies Tries to fix the above problems by the use of heuristics and are divided into two versions one hard and one soft The hard long term strategy keeps an expected minimum decrease the so called long term memory of the method If the maximum im provement of a step is smaller than this it is considered useless and t is increased instead The soft version also tries to avoid small t s but rather than forcing t to increase it just inhibits t to decrease whenever the maximum improvement is smaller than t constant Shall not be underestimated It is simple to implement and avoids the problem with too small ts but it has on the other hand no dynamic to adapt to different types and scales of problems More information can be found in Frangioni 3 pp 63 70 4 2 Performance results At the end of each example in Section 5 3 the solution and the number of iterations it took to reach convergence can be found These results can be compared with the results from the subgradient implementation described in Section 5 5 This will give some proof of the bundle implementation s s
21. hy t is called the trust region parameter The larger t the larger area around is trusted in the model 90 small t Uk Uk Figure 4 A small t gt a small step A large t gt a large step So when do we trust the model or not A simple answer is yes when a SS has occurred and no when a NS has occurred This give us the very simple rule enlarge t if a SS and reduce t if a NS Here follows the most common strategies used to update t Heuristics Here follows two of the most common heuristics used in bundle implemen tations both proposed by Kiwiel 4 The first rule compares the expected decrease with the actual decrease i e it measures how good the model is new TICO 2 9 A0 tnew is reduced when A lt 0 50 tnew is enlarged when A90 gt 0 50 there A9 6 u 0 uK 41 Observe that this is the same test as 16 the SS test with set to 0 5 The second rule is based on checking if the tentative new step is of descent 15 or ascent 2 A0 g ux 1 _ 20 2A0 g ug4 a AG tnew is reduced when g u uk 1 gt 0 tnew t tnew is enlarged when g ti ux 41 lt 0 Combinations of these two heuristics build up the most strategies used in bundle algorithms today But there are some known problems with these strategies They lack any long term memory Only the outcome of the last iterations are used to update t Only r
22. icontinuity guarantees the existence of a minimum point as soon as 0 increases at infinity 2 3 Primal dual relations How good is the dual solution in terms of the primal We know that 0 u bounds f x from above but when can we expect to obtain an optimal primal solution when the dual problem 5 is solved 2 3 1 Everett s Theorem Suppose we have maximized L u over x for a certain u and thus obtained an optimal u together with its constraint value c u R set gu c u and take an arbitrary x X such that c g By definition L x u lt L au u which can be written f x lt f tu tu u e a f au Theorem 2 Everett With the above notation u solves the following per tubation of the primal problem 2 max f z st TEXCR c x Yu If gu is small in the above pertubed problem x can be viewed as approxi mately optimal If gu 0 we are done x is the optimal solution for the primal problem 2 2 3 2 Differential study of the dual function In this section we will study the solutions of the dual function 4 We introduce the notation X u x X L x u 0 u G u g c a x E X u 7 for all solutions x for a certain u and its image through V L Because a subdifferential is always a closed convex set we deduce that the closed convex hull of G u is entirely contained in 06 u The question is if they are the same
23. ils of how to create an oracle are explained in the provided examples see Section 5 3 The easiest way to create your own oracle is to modify one of these examples Last put the two files Oracle h and Oracle C in the bundle directory 5 2 2 Parameters The file Parameters contains initial settings and can be left untouched The settings control the bundle size the t strategy the stopping conditions and some other features A full description of each parameter can be found in the file bundle h row 227 5 2 3 Main C The Main C file initiates and creates all parts of the program It contains three sections Create Creates the oracle and the solver object The x matrix contain ing all the primal solutions is also created Notice that x can at most have two indexes If you want more indexes you have to add this to the code But observe that the memory use increases exponentially with the number of indexes Solve Here is the optimal value the dual solution and the primal solution printed The primal solution is calculated as a convex combination of the saved solutions in the x matrix as in 8 The multipliers a s are given by the solver object More about this can be found in the file bundle h row 600 Delete All the created objects are destroyed 5 2 4 Makefile The file Makefile contains all the compiler directives and links all the files together with the ILOG CPLEX environment The most important rows are the
24. iterations Then we may conclude that it is unlikely that it will be part of the final model 6 that defines an optimal solution 5 The last piece to study is the stopping test Because ux 1 minimizes 0 instead of and because need not to be smaller than 0 neither nor 6 are good stop indicators But ux 41 minimizes the model 14 which is a convex function and therefore S A 1 OE 00 uK 1 O00 ux 1 SS u gt uky1 U t for some y 06 uK41 which reveals that g W Ux 1 t and can be calculated Thus we can write for all u R A u gt O u gt O uxs1 9 u ura 17 a gt 0 4 0 u ukx 1 with the conclusion that we can stop when both 6 and g are small How effective these methods are depend crucially on how t is dynamically updated which leads to t strategies studied in Chapter 4 1 More information can be found in Lemar chal 1 pp 28 30 3 6 Method connection scheme Here follows a connection scheme that reveals how the methods described are related to each other Observe how the parameter t can be used to move between them All methods except subgradient methods bundle size 1 are variations of bundle methods because they store data about previous iterations But in the literature it is often only methods based on the two subproblems 11 and 15 that are referred to as bundle methods Also notice that 11 and 15 are dual to each other and therefore in re
25. ive is initiated to an instance of IloMaximize with the Lagrangian as expression Finally we add the objective and the constraints of x to the model an instance of IloModel This model is then extracted to a solver object cplex an instance of IloCplex used when ever the problem is a linear program LP Whenever we update the model this solver object will be notified and updated as well The rest of the constructor is used to set a group of variables used by Main C and the solver object These are u size uSize m x dimensions and size xDims 1 xDim1Size n Add all unconstrained multipliers u and with an InINF ucVector new unsigned int 1 ucVector 0 InINF Is the primal problem of min type isMinProblem false 21 The most of them are obvious but the ucVector will need some explanations The bundle algorithm must know which u s are constrained and which are not All w s are constrained to R as default i e the primal problem 2 has in equality constraints This is true for this example so all we have to do is to put the stop value InInf in the first position of the ucVector If the problem has equality constraints we have to add the corresponding u s to the ucVector Example 2 will show an example of this Now we are ready to study the code solving the model This is done in the Solve function which is called every time the algorithm needs a new piece to enlarge the bundle The first thing we hav
26. multipliers ax such that L tz u 0 u for each k and LS akc zk 0 k Remark 2 Suppose the primal problem 2 has inequality constraints as in Re mark 1 Then the property 0 06 u is not a necessary condition for optimality It is instead enough and also necessary that the complementary conditions ulgu 0 u gt 0 gy gt O0 where gu 06 u are fulfilled which means that either u or g for each j has to be 0 To see why read Appendix B Proposition 1 changes in the corresponding way 2 3 3 The primal solution With the use of the x s and a s from Propostion 1 we build the primal point Crs y OL Tr 8 k Under which conditions does x solve the primal problem 2 We will give the result in three steps with increasing requirements If X is a convex set then x X Use Everett s Theorem 2 to see to what extent x is approximately optimal If in addition c is linear then cla 5 agelzk 0 k and thus 2 is feasible Because L x u f x u e x f 2 the duality gap 0 u f x tells us to what extent x is approximately optimal in 2 Remark If the primal problem 2 has inequality con straints it is enough if c is convex If in addition f is concave and thus L is concave then f a L a u gt S anL en u O u k Combine this with weak duality and x is optimal Remark 3 Observe that the algorithm solving the dual problem has n
27. n ux 1 1k 1 to minr ur eR 12 r gt 0 us 9 u u forall i The next iterate is ux 1 with its model value ux 1 rk41 A new call to the oracle with uxy1 gives a new piece O ux 1 g ux 1 to the bundle improving raising the model and the process is repeated The nominal decrease defined as 6 O uK rK41 O ux O ux41 gt 0 13 is the value that would be reduced in the true 0 if the model were accurate It can also be used to construct bounds for the optimal value 0 ux gt min 0 gt O ux 0 Tk 1 and can therefore be used as a stopping test the algorithm can be stopped as soon as 6 is small Observe that the sequence 0 ux is chaotic while the sequence rx is nondecreasing and bounded from above by any 6 u We can therefore say that rx min 8 The main problem with these methods is that an artificial box constraint has to be added to 12 at least in the initial iterations to ensure a finite solu tion Think of K 1 we get rg oo and up at infinity The cutting plane algorithm is inherently unstable and it is very difficult to find tight enough box constraints that do not cut away u More information can be found in Lemar chal 1 pp 24 26 11 3 5 Stabilized cutting plane methods To get rid of the unstable properties in 12 we introduce a stability center a a point we want ux to be near We choose a quadratic stabilizing term u
28. nfo html 33
29. o informa tion about the above requirements If X is not a convex set X will be enlarged to its conver hull It is then up to the user of Lagrangian duality with the help of the list above to analyze the optimality of the solution Remark 4 IfX is a convex set the equality constraints are linear the inequality constraints are conver and f is concave then x from 8 with a from Proposition 1 is optimal 3 Dual Algorithms In this section we will forget about the primal problem 2 and focus on methods solving the dual problem 5 First we will restate Assumption 1 with more details about the problem we have to solve Assumption 2 We are given a conver function 0 to be minimized It is defined on the whole of R The only information available on 6 is an oracle which returns 0 u and one subgradient glu for any given u R No control is possible on g u when 00 u is not the singleton VO u Remark 5 If the primal problem 2 has inequality constraints R changes to R in the assumption above and in all theory that follows in this section We will not remark on this further The third assumption above means that it is impossible to select a particular maximizer of the Lagrangian in case there are several Also notice that the assumptions above do not require the problems to be Lagrangian duals any problem that fulfills them can be solved with these algorithms
30. of linear functions These linear functions underestimate the real function but are tight in the three points u1 uz and uz Each linear function is given by the oracle as a function value 0 u and a subgradient g u for any u U 2 Lagrangian Relaxation In this section we will describe Lagrangian relaxation It is a technique to find upper bounds on an arbitrary maximization problem We will start with a rather general problem and as the theory reveals itself requirements will be added This is because we mostly are interested in problems that can be solved to global optimality with the methods we describe The theory is to a high degree based on Lemar chal 1 2 1 General technique Consider an optimization problem max f z st 2 X CR 2 cji x 0 j E 1 m where f R gt R and c R R hereafter called the primal problem We introduce the Lagrangian L x u f x gt ucila f x u cla x u E XxR 3 a function of the primal variable x and of the dual variable u R The last equality introduce the notation c for the m vector of constraint values The Lagrangian has now replaced each constraint c x 0 by a linear price to be paid or received according to the sign of uj Maximizing 3 over X is therefore closely related to solving 2 The dual function associated with 2 and 3 is the function of u defined by O u max L x u ueR 4 and the d
31. on in 5 iterations f 4 30769 x 2 15385 0 923077 1 69231 u 0 0769231 3 15385 24 5 4 Tutorials Three of the provided examples are tutorials solving two groups of problems Set Covering Problems SCP and Generalized Assignment Problems GAP For GAP we have created two different lagrangian relaxations with different performances as result To run these examples you have to follow the same list as before with some modifications Step 0 will of course not be necessary In step 2 you will not only copy the oracle files but also some data file You find other data files at 8 In step 6 will the run command change to bundle lt data file gt 5 4 1 Set Covering Problem SCP The code to this example can be found in lt Examples directory gt SCP The problem Primal problem 2 n min gt C5T j 1 s b Naya 21 for 1 m j 1 z 0 1 for j 1 n Lagrangian 3 L x u Y cia Su 1 X ait 0 1 weR i 1 j l j l Dual function 4 0 u min L a u ue R st x 0 1 The format of the data files is 1 Number of rows m number of columns n 2 The cost of each column c j j 1 n 3 For each row i i 1 m the number of columns which cover row i followed by a list of the columns which cover row 2 Results scp41 txt Optimal solution in 135 iterations f 429 small txt Optimal solution in 399 iterations f 1 6575e 06 Inc
32. oof can be found in 2 pp 295 296 However this simplicity has its price in speed of convergence and the method can only give an approximation of the optimal value More information can be found in Lemar chal 1 pp 22 24 Frangioni 3 p 5 3 2 Conjugated subgradient methods The main theoretical problem with subgradient methods is that subgradients are not guaranteed to be directions of descent Conjugated subgradient methods try to detect this by making an exact line search LS in the current direction If this direction is not of descent LS will just return the current point u During such a LS a new subgradient in a neighbourhood of u can be generated this possibility suggests making some form of inner approximation of the subdiffer ential where the bundle is meant to contain some subgradients of 0 at the current point A new direction can then be found by minimizing the vector norm over the convex hull of the known subgradients at a min Y nyil p O min lgl g conv gi 4 B 9 i B where y ies yi 1 gt 0 is the unit simplex in the G dim space and the new direction is d gt a HiPi The result of the iteration is thus either a significantly better point u or a new subgradient g which is not guaranteed to be a subgradient in u but to yield a new direction in the next iteration since it is obtained close to u it still gives relevant information about 0 around uw
33. ound in the file Parameters Results File File size Iterations Optimal value heuristic soft hard scp45 txt 20 kB 64 64 123 512 scp51 txt 44 kB 173 173 257 251 225 scp61 txt 44 kB 225 225 228 133 14 scpal txt 95 kB 772 759 437 246 836 scpcl txt 167 kB 723 624 317 223 801 As can be seen from these results the long term strategies especially the hard one seem to be best on large problems and the pure heuristic seems to be best on the smaller problems 32 References 1 C Lemar chal Lagrangian Relaxation Computational Combinatorial Op oOo N QO Oo timization M Juenger and D Naddef eds Lecture Notes in Computer Science 2241 115 160 Springer Verlag 2001 T Larsson M Patriksson A B Str mberg Ergodic primal convergence in dual subgradient schemes for conver programming Math Program 86 283 312 Springer Verlag 1999 A Frangioni Dual Ascent Methods and Multicommodity Flow Problems Ph D Thesis Part 1 Bundle Methods 1 70 1997 K C Kiwiel Proximity Control in Bundle Methods for Convex Nondiffer entiable Optimization Math Program 46 105 122 1990 ILOG CPLEX 8 1 Getting Started 2002 ILOG Concert Technology 1 3 User s Manual 2002 ILOG Concert Technology 1 3 Reference Manual 2002 J E Beasley OR Library distributing test problems by electronic mail Journal of the Operational Research Society 41 11 1069 1072 1990 website http mscmga ms ic ac uk i
34. rease the maximum value for t parameter in the file Parameters and see how this effects the number of iterations More results can be found in Appendix D 25 5 4 2 Generalized Assignment Problem GAP The code to this example can be found in lt Examples directory gt GAP 1 2 The problem Primal problem 2 min 5 y Cig Vij i 1 j 1 m s b Sg for j 1 n i 1 n eae for 1a j l xi E 0 1 for i 1 m j 1 n Dual function 4 GAP1 i 1 j 1 j 1 O u min L x u min 5 Cig Vig Sou Tij 1 j i l n s t ayy lt ri for i 1 m j 1 xi E 0 1 for i 1 m j 1 n Dual function 4 GAP2 m n m n O u min L x u min gt Cig Lig gt Uj gt QijZij fi i 1 j 1 i 1 j 1 m s t Y ty 1 for j 1 n i 1 vig 0 1 for 1 mj 1 n As can be seen the two relaxations results in two different problems the first one can be recognized as the knap sack problem and the second one as the simple assignment problem Which one is to prefer One usually say that you shall not relax too much of the original problem and in that case the first relaxation is the preferred one Can this be seen in the result The format of the data files is 1 Number of agents m number of jobs n 2 For each agent i i 1 m in turn cost of allocating job j to agent i j 1 n 26 3 For each agent i i 1 m in turn re
35. se System SYSTEM ultrasparc_5_5 0 LIBFORMAT static_pic_mt CPLEX directories CPLEXDIR usr site pkg cplex 8 1 cplex81 CONCERTDIR usr site pkg cplex 8 1 concert13 The first two describe the system and the last two contain the paths to ILOG CPLEX and ILOG Concert Technology If the system is updated or a new version of CPLEX is installed these settings have to be changed If so try to find the correct Makefile for all examples CPLEX sends along each system and version The current ones can be found in lt cplex81 directory gt examples lt machine gt lt libformat gt 19 5 3 Examples In this section we will present some examples showing how to create oracles with ILOG Concert Technology The first example will be explained in detail while only the new things will be highlighted in the second one All the files are found in the Example directory unpacked in Section 5 1 The examples will also be a good starting point when creating your own oracles It s good to have access to the CPLEX manuals 6 and 7 while studying the files Any code starting with Ilo is either a class or function in the Concert Technology library 5 3 1 Example 1 The code to this example can be found in lt Examples directory gt EX1 The problem Primal problem 2 max 121 222 s t ei a 4x9 lt 8 0 lt xr lt 4 Lagrangian 3 L a u x 2 2 u 11 4z2 8 0 lt x lt 4 ueR Dual function 4 0 u m
36. source consumed in allocating job j to agent i j 1 n resource capacity of agent i i 1 m Results GAP1 gap11 txt Near optimal solution in 48 iterations f 337 GAP2 gap11 txt Solution in 18 iterations f 343 587 The optimal solution is 336 The first problem is harder to solve but gives a better result The answer is thus yes to our question above 5 5 Subgradient implementation To be able to make comparisons between the bundle implementation and a sim ple method we have implemented the subgradient method described in Section 3 1 Follow these steps to run it 1 Unpack the file Subgradient tar into a directory of your choice tar xvf Subgradient tar 2 Create an oracle Oracle h and Oracle C or copy one of the five pro vided examples Put the files in the Subgradient directory See Section 5 2 1 3 Optionally modify Main C to change the step size rule and printouts The rule used is the one described in Section 3 1 5 Compile the program with the command make 6 Run the program with the command sg As you can see you use it in the same way as the bundle implementation with the same oracles as problem files Results after 1000 iterations EX1 f 6 04031 x 3 96 1 u 0 506076 EX2 f 4 25305 aw 2 11 0 906 1 694 u 0 084553 3 15582 SCP scp41 txt f 427 264 SCP small txt f 52 3913 GAP1 gapii txt f 357 375 GAP2 gap11 txt f 343 742
37. ual problem is then to solve in 0 u 5 win 0 u 5 Lagrangian relaxation accepts almost any data f X c and can be used in many situations One example is when we have both linear and nonlinear con straints some of them covering all variables making it impossible to separate The only crucial assumption is the following it will be in force throughout Assumption 1 Solving 2 is difficult while solving 4 is easy An oracle is available which solves 4 for given u R Lagrangian relaxation then takes advantage of this assumption and aims at finding an appropriate u To explain why this is the same as solving the dual problem 5 observe the following by the definition of 6 0 u gt f x for all z feasible in 2 and all u R simply because u c w 0 This relation is known as weak duality which gives us two reasons for studying the dual problem O u bounds f x from above Minimizing 0 amounts to finding the best possible such bound For the dual function 4 to be any good it must produce some argmax Zu Of L u which also solves the primal problem 2 This u must be feasible in 2 and therefore satisfy 0 u f x From weak duality we see that the only chance for our u is to minimize 6 Remark 1 Suppose the primal problem 2 has inequality constraints max f a st 2 X CR cj x lt 0 j 1 m Then the dual variable u must be nonnegative an
38. uperiority over simpler methods More results about the performance of the bundle algorithm can be found in Frangioni 3 pp 46 50 pp 69 70 16 5 Manual The bundle implementation is built around professor Antonio Frangioni s 3 code developed during his Ph D research It s written in C and is highly advanced with many settings not covered in this report The code is built to work on a UNIX system and demands that ILOG CPLEX version 8 or newer is installed Here follows a list of steps necessary to solve a problem with this bundle implementation 0 Relax the problem you want to solve See Section 2 1 1 Do some preparations See Section 5 1 2 Create an oracle Oracle h and Oracle C or copy one of the five pro vided examples Put the files in the Bundle directory See Section 5 2 1 Optional modify initial parameters See Section 5 2 2 Optional modify Main C to change the printouts See Section 5 2 3 Compile the program with the command make Run the program with the command bundle 5 1 Preparations You will need the following files Bundle tar and Examples tar Start by un packing them into a directory of your choice tar xvf Bundle tar tar xvf Examples tar This will create two directories containing the bundle implementation and the examples explained in Section 5 3 The CPLEX environment needs to be initialized once before its use setenv ILOG_LICENSE_FILE usr site pkg cplex 8 1 ilm ac
Download Pdf Manuals
Related Search