# Applications of Linear/Integer Programming

### IBM OSL Studies

The following excerpts are from the EKKNEWS newsletter describing various uses of IBMs Optimization Subroutine Library (OSL).

## Emergency Electric Power Restoration Using OSL

by Bhaba R. Sarker and Srinivas E. Mahankali

The present work is a part of the ongoing research project at the Department of Industrial & Manufacturing Systems Engineering at Louisiana State University, sponsored by Louisiana Education Quality Support Funds (LEQSF).

The problem of assigning electric power repair depots to various locations in a damaged area assumes critical importance during emergency situations such as natural disasters. This issue has to be addressed taking into account the damage levels, priorities, and the consequent demand for various resources among different cells in the area, and also the capacity considerations for the depots.

This is a typical depot allocation and vehicle routing problem with several other prevailing factors such as repair time and time windows for repair, crew and resource scheduling with resources' precedence contraints and demands for same resource from multiple locations (cells in a grid mapping system). Considerable amount of literature exists in the area of location research, it may be mentioned that the problem of simultaneous location of depots with varied capacities, and their assignment to customers with different demands has not been dealt with. So all these constraints make the problem immediately complex and the mathematical structure tends to be forbidden.

For the purpose of power restoration, a region may be divided into a number of cells. Each cell may be differently affected by the disaster, and the damage levels may be assessed to ascertain the amount of resources required for power restoration in each cell. In a similar vein, the available depots may also differ in terms of their capacities with respect to various resources. The problem, then, is to determine the location of the available depots and also, to ascertain the amount of resource transported from different depots to various cells, so as to minimize the total transportation cost. It is implicitly assumed that this objective also enables restoration of power in minimal time.

As mentioned earlier, the demand arising at various cells is dependent on the level of damage and priorities which are different for each cell. This, coupled with the fact that all the depots may not be identical in terms of their capacity, leads to formulation of the depot allocation problem as a mixed integer quadratic programming problem, which has been developed during the research work.

The mixed integer quadratic Programming model has a term involving the product of an integer (0-1) variable and a continuous variable. The constraints essentially reflect the demand requirement at different cells, capacity of depots for various resources, and the requirement that each depot may be located at only one cell and each cell may in turn contain only one depot. The objective function, as mentioned earlier, is a cost function which is to be minimized.

The above problem has been reformulated as an equivalent mixed integer linear programming problem. In the latter, the objective function contained only the real variables. The modification has been effected by introducing a new variable to replace the product of the integer and real variables in the previous formulation, and adding another set of constraints which impose the equivalence of these two terms, by the well known "big M" method. To have an idea of the problem size, for a problem involving with m resources, N depots, and n cells in the region, the total number of variables is Nn(1+ m + mn), of which there are Nn number of integer variables. The number of constraints is given by (n + N + mn + 2Nmn).

The datasets were randomly generated and the solution is obtained using the mixed integer programming solver in IBM's Optimization Subroutine Library (OSL) and GAMS software package on a RISC 6000 (Model 590) workstation. The optimization problem is solved for various dimensions of the problem in terms of the number of cells in the region, the number of depots considered, and the number of resources used.

Our driver program is a relatively straight-forward one. It reads in model data in the MPS format from a file (EKKMPS). Then, the read-in model is scaled (EKKSCAL) to avoid any numerical instabilities that may be encountered during the solve process. This is followed by a "presolve" phase (EKKPRSL) of the resulting scaled matrix which attempts to remove any redundancies by eliminating rows and fixing variables. This phase results in transforming the model into a new set of rows and columns. The crash routine (EKKCRSH) is then invoked to identify a starting basis for solving the LP relaxation of the MIP model. OSL's primal simplex solver (EKKSSLV) is then called to obtain an optimal solution to the relaxed problem. If the resulting solution does not satisfy the integrality requirements, then OSL's preprocessor EKKMPRE is called to further investigate the underlying polytope and generate feasible solutions quickly. It is after this phase that the branch-bound-solver in OSL (EKKMSLV) is invoked to obtain an optimal solution to the problem. After obtaining an optimal integer solution, the model is "postsolved" (EKKPSSL) to transform the model data back to the original variables and rows.

While the OSL functions called from our driver program remained the same for all our test problems, our experimentation involved manipulating various control variables to influence the solution procedure. Our computational experience shows that OSL entails immense flexibility in terms of providing the user with various user-desired options, and the performance of the software was very satisfactory in terms of the computational run times (see table below).

` Problem Number of Number of Number of CPU Time Rows Columns 0-1 Variables (in seconds)  a1 327 312 12 2.45 b1 405 396 18 8.30 a2 669 648 18              15.06 b2             683           672          24              28.52 a3            1131          1104          24              61.08 b3            1033          1020          30             104.60 b4            1455          1440          36             259.90 a4            1713          1680          30             200.20 b5            1949          1932          42             538.16 a5            2415          2376          36             556.01 b6            2515          2496          48             1028.54 a6            3237          3192          42             1119.29  `

#### References

1. Brooke, A, Kendrick, D. and A. Meeraus, 1988. GAMS: A User's Guide, The Scientific Press, Palo Alto.
2. Jayaraman, V.J. and R. Srivastava, 1995. A service logistics model for simultaneous siting of facilities and multiple levels of equipment. Computers and Operations Research 22, 191-204.
3. Nemhauser, G.L. and L.A. Wolsey, 1988. Integer and Combinatorial Optimization, John Wiley.
4. Schrage, L, 1986. Linear, Integer and Quadratic Programming with LINDO, The Scientific Press, Palo Alto.
5. Williams, H.P., 1978. Model Building in Mathematical Programming, John Wiley.

## OSL Applied in the Strategic and Tactical Management of Oil Spill Clean Up

by Wilbert E. Wilhelm, Anand V. Srinivasa and Gyana R. Parija

The environmental impact of oil spills can be dramatic as evidenced by just one recent spill. The Exxon Valdez spill occurred in 1989 and Exxon has paid over \$4 billion for clean up and ecological damage since then (Hann, 1992). Recently, a federal jury found Exxon liable for an additional \$5 billion to compensate fishermen and others who claimed economic losses from environmental damage (U.S. News and World Report, 1994). Over 60% of the crude oil imported into the U.S. travels through Texas Gulf region (Rainey, 1992), and during 1978-1992 tanker spills surpassed those in all other U.S. regions (Schneider, 1992). Spills are often associated with tanker traffic but can originate from numerous sources including oil well blow outs, aging pipelines, and underwater seepage (Geyer, 1981). Texas led the nation in oil spills in 1990 (Rainey, 1992) and remains especially vulnerable.

Effective clean up can be enabled by a hierarchial approach, consisting of strategic and tactical levels. Strategic decisions deal with investment and contingency planning strategies. The Oil Pollution Act of 1990 (OPA-90) requires an Area Contingency Plan to enable an effective and timely response to any oil spill likely to occur (U.S. Coast Guard, 1992). Tactical decisions, which are made after a spill has occurred, prescribe response systems necessary to meet the time-phased clean up requirements stipulated by law.

Currently, four methods are used to clean up oil from the sea: mechanical containment, chemicals (e.g., dispersants), bioremediation, and burning. A number of factors determine the efficacy of a method, including the amount and type of oil; the ecological, oceanographic, and climatologic conditions; and the season. An effective method, or combination of methods, must be prescribed to deal with the characteristics of a given spill. Clean up capability is provided by sy stems composed of components (e.g., equipments, materials, and personnel) at staging areas for deployment on scene.

Our research (Srinivasa and Wilhelm, 1994; Wilhelm et al, 1994; Wilhelm and Srinivasa, 1994) formulated the strategic and tactical decision problems as integer programs and developed effective solution methods for both. The models incorporate practical features, accounting, for example, for component availability and degradation of system capabilities over time as the oil weathers in the sea. Both successfully solved large-scale test problems that represent application to the Galveston Bay Area (Wilhelm et al, 1993).

The strategic, area-wide contingency planning model (Wilhelm et al, 1994) accepts as inputs the set of risk points and the likely spill scenarios and response requirements for each, the sites of existing storage locations and the inventory of components at each, and potential sites for new storage locations. It prescribes a minimum total cost plan to build new storage locations and/or expand existing ones and to purchase new components and pre-position them, and a contingency plan which determines what systems should be composed to enable an effective time-phased response for each likely spill scenario. A family of heuristics based on linear programming (LP) relaxation (Wilhelm et al, 1994) resolve the strategic model, prescribing an area-wide contingency plan. Computational tests indicate that four heuristics are quite effective, prescribing solutions for each of 10 test cases within 1.41 percent of optimum and within a few minutes runtime.

The objective of the strategic model is to minimize the present worth of total cost, including the fixed investment for opening new warehouses or expanding existing ones, the variable cost for additional storage space at new and expanded warehouses, the cost of purchasing and maintaining new components, the expense of relocating existing components, and the cost of composing and operating response systems. Five types of decision variables are involved: Z variables determine if a new warehouse should be opened at each potential new location or if an exiting one should be expanded; F variables prescribe the storage space for each new and expanded warehouse; R variables detail the relocation of existing components; X variables specify what new components should be purchased and where each should be stored; and Y variables give the contingency plan, specifying the number of response systems and their deployment times for each likely spill scenario. The model is composed of four types of constraints in addition to those that specify variables: F variables are real; Z variables are binary; and R, X, and Y variables are non-negative integers. The first type of constraint assures that capacity is provide d at a location only if a new warehouse is opened or an existing one is expanded there. The second type sets the capacity of each new or expanded warehouse large enough to accommodate the components stored there. The third makes sure that response systems use no more than the number of components made available either from existing inventories or from purchases. Finally, the fourth type of constraint assures that response systems will provide required, time-phased response to spills of all likely scenarios. The large-scale test problems that represent application to the Galveston Bay Area entail some 1869 variables and 3264 constraints. The solution procedure consists of various pre-processing methods, two heuristics and an optimizing approach. One of the pre-processing methods exploits Generalized Upper Bound (GUB) constraints in the model to tighten the formulation. This is achieved by iteratively solving the Linear Programs (LP) obtained by replacing the original objective function with the LHS of each GUB constraint. The LPs were solved using OSL's LP solver EKKSSLV. Other methods to tighten the formulation like obtaining bounds on the number of response systems necessary also relied on EKKSSLV. Initial efforts to obtain feasible solutions using a straight forward call of EKKMSLV were discouraging. For some test cases, we failed to obtain a feasible solution even after 10 hours of run time. Hence, heuristics were devised with the goal of obtaining 'good' feasible solutions quickly. The first heuristic is totally based on LP relaxation and obtains feasible solutions by judiciously setting lower bounds on certain variables. The choice of which variables to set relied on a property, 'response system splitting.' Due to the degradation of response systems' capabilities over time once a system is deployed, the LP solution typically prescribes splitting whereby the same response system is deployed at different time points, each time at a fractional value, but such that all fractional values sum exactly to an integer. Splitting means that the heuristic can set the appropriate lower bounds on variables representing the deployment of a system at the time points involved without violating feasibility vis-a-vis resource availability, because the new bounds merely sum to the total response prescribed by the LP relaxation. The second heuristic combines LP relaxation with partial enumeration to obtain a feasible solution. The procedure completes one iteration of Heuristic I to tighten variable lower bounds and then invokes OSL's integer solver EKKMSLV. Both the heuristics have been very successful in obtaining 'good' feasible solutions quickly. This example shows how efficient pre-processing methods for general integer variables can be devised by exploiting structures in the problem and using the solution capabilities of OSL. Furthermore, the performance of EKKMSLV can be dramatically improved by exploiting an in-depth understanding of the problem being addressed. Inputs to the tactical decision problem (Srinivasa and Wilhelm, 1994; Wilhelm and Srinivasa, 1994) include the location of the spill and its characteristics, the time-phased response required, the inventory of components at each storage location, and the sites and capacities of potential staging areas. The model prescribes a minimum time response, guarding against environmental damage by dispatching components from their storage locations to staging areas where appropriate systems can be composed for timely deployment on scene. Our solution approach includes a unique column generation method to compose systems and several pre-processing methods to facilitate calculations by bounding decision variables. A heuristic (Wilhelm and Srinivasa, 1994) based on LP relaxation found solutions for 12 test cases within 2.9 percent of optimum in less than 1.1 seconds (on average). A new optimizing method (Srinivasa and Wilhelm, 1994), which based pre-processing methods and generation of families of facets on an aggregated model, solved each of 6 cases in less than 15 minutes and each of 12 test cases in less than 1.75 hours.

The objective of the tactical decision problem is to minimize total response time, reflecting the consensus that, once a spill has occurred, it is most important to respond as quickly as possible to ameliorate environmental damage. Decision variables prescribe the number of response systems of each type that must be deployed at each time point. The model consists of three types of constraints in addition to the non-negative integer requirements on the decision variables. The first type of constraint assures that response systems provide the required, time-phased response. The second type assures that deployment of response systems is limited by the components available, and the third type imposes restrictions on composing response systems resulting from the limited capacities of staging areas. The large-scale test problems that represent application to the Galveston Bay Area entail 180 constraints and the column generation method generated some 350 columns, specifying the components that compose each system, the storage location of each component, and the staging area where each system is composed. The optimizing approach is based on pre-processing methods and the generation of families of facets of the polytope associated with an aggregated model of the TDP. The aggregated model which is also a general integer program but of substantially smaller size (the model has about 75-85% fewer variables), is solved using EKKMSLV. However, instead of assuring the integrality of decision variables using EKKMPS, we imposed integrality by using EKKIMDL. This was done to take advantage of the flexibility of EKKIMDL in defining SOSs by placing each general integer variable in its own set. This allowed us to define priorities for each set (i.e., each general integer variable) based on the clean up capability and time availability of the response system represented by that variable. Such a scheme for prioritizing the decision variables provided for a convenient branching and resulted in a significant improvement in the performance of EKKMSLV.

These methods can prescribe area-wide contingency plans and tactical responses for actual cases in reasonable time. This capability should enable dramatically improved oil spill clean up response.

This material is based on work supported by the Texas Advanced Technology Program under Grant Number 999903-282.

References.

1. Geyer, R.A., Naturally Occurring Hydrocarbon Seeps in the Gulf of Mexico and the Caribbean Sea, Dept. of Oceanography, College of Geosciences, Texas A&M Univ., 1981.
2. Hann, R. W., Course Manual: Oil Spill Prevention and Control for the 90's, Environmental Engineering Program, Civil Engg. Dept., Texas A&M Univ., College Station, Texas, 1992.
3. Rainey, G., "The Risk of Oil Spills from the Transportation of Petroleum in the Gulf of Mexico," Proceedings: Seminar on Physical Recovery of Spills, hosted by the Eighth Coast Guard District, the State of Texas, the Port of Corpus Christi, and the Corpus Christi State Univ., Corpus Christi State Univ., Corpus Christi, July 27-8, 1992, 131-142.
4. Schneider, M., Houston Chronicle, August 24, 1992.
5. Srinivasa, A.V., and W.E. Wilhelm, "An Optimizing Procedure for Tactical Response in Oil Spill Clean Up Operations," Working Paper, Dept. of Ind. Engg., Texas A&M Univ., College Station, TX, 1994.
6. United States Coast Guard, Commandant Notice 16471, 1992.
7. U.S. News and World Report, "Exxon's New Valdez Bill: \$454 per Gallon," Sept. 26, 1994, p. 29.
8. Wilhelm, W.E., Geyer, R.A., and S. Joo., "A Manual for Formulating Data for the Strategic Oil Spill Clean Up Model," Working Paper, Dept. of Ind. Engg., Texas A&M Univ., College Station, TX, 1993.
9. Wilhelm, W.E., Parija, G.R., and A.V. Srinivasa, "A Strategic, Area-Wide Contingency Planning Model of Oil Spill Clean Up Operations," Working Paper, Dept. of Ind. Engg., Texas A&M Univ., College Station, TX, 1994.
10. Wilhelm, W.E., and A.V. Srinivasa, "Tactical Response in Oil Spill Clean Up Operations," Working Paper, Dept. of Ind. Engg., Texas A&M Univ., College Station, TX, 1994.

## Stochastic Programming at IBM's T.J. Watson Research Center

by Alan J. King

Linear programming under uncertainty, or more generally, stochastic programming, is a field reborn in the new computing frontier of powerful processors networked together. Software to test approaches to modeling and solving stochastic programs is under development at IBM's Thomas J. Watson Research Center, and has proven of crucial value in early applications of multiperiod risk management in the finance industry. As Roger Wets said in his address in honor of George Dantzig at the recent Mathematical Programming Symposium in Ann Arbor, from the earliest days of linear programming Dantzig held that stochastic programming was "the real problem". Real data is always uncertain. Real decisions always have uncertain outcomes. The real problem always involves risk. In many real world situations, there are various causes of uncertainty - interest rates may change; demand may increase or decrease; tax laws may change; a merger may or may not occur. A "scenario" is the term used to refer to the situation we get when each of the causes of uncertainty is fixed to a single value. Clearly some scenarios are more likely than others! An example from financial portfolio management illustrates the issues. For a given scenario of financial market returns, one could find a portfolio composition that maximizes shareholders equity subject to regulatory constraints. But imagine yourself as the CEO of a large insurance company analyzing 100 different scenarios of possible financial market returns. Each one comes with a different portfolio recommendation. In one scenario, it may be best to invest in stocks. In another, to invest in bonds. In yet another, it may be best to get out of town! But you have to make a decision now, while everything is still uncertain. The CEO in our example would be looking for solutions that could be implemented today and that do the best job with respect to all the risks faced by the company. The problem, of course, is that scenario analysis gives solutions that would be optimal under that scenario, only. Scenario solutions are prescient: they know exactly what is going to happen tomorrow and the day after that. Furthermore, the scenario solutions don't agree on what to do today. Solutions for one scenario can turn out to be disastrous under another scenario. And solutions that minimize risk, taking into account the full spectrum of possible scenarios, sometimes don't show up as optimal in any one of the scenarios. Stochastic programming is a methodology for risk management. To help the CEO analyze the risk of a decision made today, two things need to be done. First, the scenarios must be combined to give a stochastic multistage distribution, so that solutions made in one time stage will have to hedge against the uncertainty in the next stage. Second, the risk of a decision must be formed by computing costs or penalties under each scenario. For example, when a target is missed or an expensive recourse action must be undertaken, the penalty or cost is incorporated into the risk function. Stochastic programming seeks a solution that minimizes this risk. For a decision maker facing the given distribution of uncertain scenarios and who possesses the given risk profile, such a solution would be the optimal hedge. A single instance of a stochastic program is huge -- typical sizes for practical models run into many hundreds of thousands of rows. Approximation techniques and specialized decomposition algorithms are absolutely necessary for the solution of problems of such enormous scale. Of course, in most situations one does not have a precise idea of one's risk function or the stochastic distribution of uncertain events. In a single decision situation, many different risk objectives would be selected and compared, and many scenarios would be generated and combined into different distributions. Even more than in linear programming, parametric programming and sensitivity analysis forms a crucial part of stochastic programming methodology. Only recently have computers attained the minimal power and performance necessary to support risk management at the level of detail practical for industrial application. Nevertheless, stochastic programming is at a stage now where it can make a serious impact on decision situations. Research software, dubbed SP/OSL, under development by Alan King and colleagues in the Optimization and Statistics Center at IBM's Thomas J. Watson Research Center, implements a stochastic programming interface using OSL routines as solution modules. SP/OSL supports basic stochastic programming activities, such as bundling scenarios and generating stochastic programs with linear or piecewise linear risk functions. It also incorporates decomposition tools for quick solution. The insurance company asset allocation problem discussed above is an actual example from a consulting engagement led by Alan King. The company manages more than \$20 billion in assets. A multi-year cash flow model linking portfolio transactions and earnings and underwriting losses to basic accounting measurements, like shareholders equity and internal rate of return, was built by members of IBM consulting practices under the guidance of insurance company experts. Scenarios of asset returns, inflation, and underwriting losses were generated from the company's proprietary forecasting models. The scenarios were bundled by SP/OSL routines, and risks modeled as penalizations above or below company-selected targets. The stochastic programs solved by SP/OSL for the insurance company include many more scenarios than the company's scenario analysis procedure ever did. SP/OSL solves 540-scenario stochastic programs in less than an hour. (A 540-scenario problem has about 350,000 rows.) Furthermore, the exploration of the company's risks and opportunities is far more sophisticated and comprehensive than before. The CEO is now presented with a chart summarizing the output of parametric runs comprising the solutions of over 100 stochastic programs representing different blends of risk measurements and assumptions. An early version of SP/OSL was instrumental in demonstrating the practicality of stochastic programming for the prize-winning Russell-Yasuda Kasai asset/liability management model (see reference 1). SP/OSL is now equipped to handle many types of distributions and to interface with LP modeling languages. SP/OSL resembles OSL, in that it is a library of subroutines that can perform various actions on stochastic program data and models. The decomposition method built into SP/OSL has recently been parallelized (by its designer, Steve Wright) using PVM subroutines. Several universities and research labs are field-testing SP/OSL for use in risk management applications and/or implementation of solution methods, including: University of Geneva, University of Michigan, Georgia Tech, University of California at Davis, University of Texas at Austin, and the University of Cyprus. Applications under development include: long-range energy planning, reservoir management, commercial aircraft arrivals scheduling, and production/inventory management, to name just a few. Reference 2 documents the present state of the software. Interested readers may contact Alan King (kingaj@watson.ibm.com) for further information.

References.

1. David R. Carino, et al, "The Russell-Yasuda Kasai Model: An Asset/Liability Model for a Japanese Insurance Company Using Multistage Stochastic Programming", Interfaces 24 (1994) pp. 29-49.
2. Alan J. King, Stephen E. Wright, and Robert Entriken, "SP/OSL Version 1.0: Stochastic Programming Interface Library User's Guide", (Technical Report, IBM Research Division, Yorktown Heights, New York, 1994).

## An OSL-Based Implementation of a Cutting Plane Method for Solving Assembly System Design Problems

by Gyana R. Parija and Wilbert E. Wilhelm