1. Introduction
Consider the linear programming problem
(P)
(1)
s.t.
(2)
(3)
where x and c are the n-dimensional column vectors of variables and objective coefficients, respectively, and z represents the objective function. The matrix A is an m × n matrix 
with row vectors
; b is an m-dimensional column vector; and 0 is an n-dimensional column vector of zeros. The non-polynomial simplex methods and polynomial interior-point barrier-function algorithms illustrate the two different approaches to solve problem P. There is no single best algorithm [1] . For any existing approach, there is a problem instance for which the developed method performs poorly [2] [3] . However, interior point methods do not provide efficient post-optimality analysis, so the simplex algorithm is the most frequently used approach [2] , even for sparse large scale linear programming problems where barrier methods perform extremely well. In fact, the simplex method has been called “the algorithm that runs the world” [4] ; yet it often cannot efficiently solve the large scale LPs required in many applications.
In this paper we consider both the general linear program (LP) and the special case with 
and
,
;
; and
, which is called a nonnegative linear program (NNLP). NNLPs have some useful properties that simplify their solution, and they model various practical applications such as determining an optimal driving route using global positioning data [5] and updating airline schedules [6] , for example. We propose active-set methods for LPs and NNLPs. Our approach divides the constraints of problem P into operative and inoperative constraints at each iteration. Operative constraints are those active in a current relaxed subproblem Pr, 
of P, while the inoperative ones are constraints of the problem P not active in Pr. In our active-set method we iteratively solve Pr, 
of P after adding one or more violated inoperative constraints from (2) to 
until the solution 
to Pr is a solution to P.
Active-set methods have been studied by Stone [7] , Thompson et al. [8] , Myers and Shih [9] , and Curet [10] , for example. More explicitly, Adler et al. [11] suggested a random constraint selection rule, in which a violated constraint was randomly added to the operative set. Bixby et al. [1] developed a sifting method for problems having wide and narrow structure. In effect, the sifting method is an active-set method for the dual. Zeleny [12] used a constraint selection rule called VIOL here, which added the constraint most violated at each iteration. Mitchell [13] used a multi-cut version of the VIOL for an interior point cutting plane algorithm. Corley et al. [14]
developed a cosine simplex algorithm where a single violated constraint maximizing 
was
added to the operative set. Junior and Lins [15] used a similar cosine criterion to determine an improved initial basis for the simplex algorithm.
References [6] [16] and [17] are the most directly related to the current work. The constraint optimal selection technique (COST) RAD was introduced in [6] and [16] as a constraint selection metric for NNLPs, and then generalized in [17] to GRAD for LPs. In RAD and GRAD, an initial problem P0 is formulated from P such that all variables are bounded by at least one constraint, which may be an artificial bounding constraint ![]()
for sufficiently large M. Multiple violated constraints are then added to problems Pr, ![]()
according to the constraint selection metric RAD or GRAD. The contribution here is to improve significantly the speed of RAD and GRAD by varying the number of added constraints at each Pr according to an empirically derived function estimating the effectiveness of that iteration.
The paper is organized as follows. The constraint selection metric and dynamic active-set conjunction with RAD and GRAD is given in Section 2. In Section 3, we describe the problem instances and CPLEX preprocessing settings. Section 4 contains computational experiments for both NNLPs and LPs. Then the performance of the new methods is compared to the previous ones as well as the CPLEX simplex, dual simplex, and barrier method. In Section 5, we present conclusions.
2. A Dynamic Active-Set Approach
The purpose of our dynamic active-set method is to add violated constraints to problem Pr more effectively than in [6] and [17] . In COST RAD of [6] for NNLPs we use the constraint selection metric
![]()
(4)
to order constraints from highest to lowest value of RAD as a geometric heuristic for determining the constraints most likely to be binding at optimality. Moreover, at each iteration in [6] we add violated constraints in order of decreasing RAD until the added constraints contain non-zero coefficients for all variables. In similar fashion for COST GRAD of [17] , we use the constraint selection metric
![]()
(5)
where
![]()
(6)
for a small positive constant ε and ![]()
However, in COST GRAD, violated constraints are added until every column of A has at least one positive and one negative value.
In this paper we propose a dynamic method that adds a varying number of constraints to Pr that depends on the progress made at![]()
. No equality constraints are considered here, but any equality constraints can be included in P0. An active-set function is defined to compensate for the lack of progress in ![]()
by adding more violated constraints at Pr. The algorithm stops when the solution ![]()
to Pr is the optimal solution to P. Of course, one could simply add all violated constraints at any one iteration. However, the proposed dynamic method attempts to balance the number of iterations required to eliminate all constraints violation, while still maintaining an efficient active-set method.
2.1. Dynamic Active-Set Approach for NNLP
The dynamic active-set approach developed for solving NNLPs is as follows. Constraints are initially ordered by the RAD constraint selection metric (4). To construct P0, we choose constraints from (2) in descending order of RAD until each variables xj has an ![]()
in the coefficient matrix of P0. We say the variables are covered by the constraints of the initial problem P0. P0 is then solved by the primal simplex to achieve an initial solution![]()
. Now let ![]()
be the number of constraints in the original problem P and, in general, let ![]()
be the number of constraints of problem P violated by![]()
. At each iteration r, the total number of violated constraints ![]()
is computed, and the improvement percentage is calculated by
![]()
(7)
where ![]()
represents the percent of improvement made in reducing the total number of violated constraints at iteration r. Next, with [.] denoting the greatest integer function, let
![]()
(8)
where![]()
. The value ![]()
is an upper bound on the possible number of non-operative violated constraints that can be added at active-set iteration ![]()
The actual number added is![]()
. The active-set iterations terminate when ![]()
and therefore![]()
. Equation (8) was developed from the results of com- putational experiments.
Pseudocode for the dynamic active-set NNLPs is as follows.
Step 1―Identify constraints to initially bound the problem.
1: ![]()
2: while ![]()
do
3: Let ![]()
4: if ![]()
and ![]()
then
5: ![]()
6: end if
7: ![]()
8: ![]()
9: end while
Step 2―Using the primal simplex method, obtain an optimal solution ![]()
for the initial bounded problem P0
![]()
![]()
![]()
Step 3―Perform the following iterations until an optimal solution to problem P is found.
1: ![]()
2: ![]()
3: ![]()
4: ![]()
5: ![]()
6: ![]()
7: Calculate ![]()
8: ![]()
9: ![]()
10: for (i = 0 to![]()
) ![]()
11: Solve the following Pr by the dual simplex method to obtain ![]()
12: else if (ωr = 100) then Optimized ¬ true//![]()
is an optimal solution to P.
13: end if
14: else Optimized ¬ true//![]()
is an optimal solution to P.
15: end if
16: end while
2.2. Dynamic Active-Set Approach for LP
The dynamic active-set approach for solving LP similar to the one for NNLPs. We construct P0 by choosing a number of constraints ![]()
from (2) in descending order of GRAD from (2) until no variable xj is left without at least either a positive or a negative coefficient. In addition, we include an artificial bounding constraint![]()
. If![]()
, then set![]()
. Then P0 is solved to obtain an initial solution![]()
. As in Section 2.1, it is initially assumed that all constraints are violated![]()
. Then the relative improvement percent ![]()
is calculated by (7) for Pr and![]()
. Now let
![]()
(9)
where the value ![]()
is an upper bound on the possible number of non-operative violated constraints that can be added at active-set iteration ![]()
The actual number added is![]()
. As ![]()
decreases, ![]()
increases in (9) to add more violated constraints to![]()
. The algorithm stops at 100% reduction in the number of violated constraints.
Pseudocode for the dynamic active-set for LPs is as follows.
Step 1―Identify constraints to initially bound the problem.
1: ![]()
2: while ![]()
do
3: Let ![]()
4: if ![]()
and ![]()
or ![]()
then
5: ![]()
6: end if
7: ![]()
8: Optimized ¬ false
9: end while
Step 2―Using the primal simplex method, obtain an optimal solution ![]()
for the initial bounded problem P0 given by
![]()
subject to ![]()
![]()
![]()
Step 3―Perform the following iterations until an optimal solution to problem P is found.
1: ![]()
2: ![]()
3: ![]()
4: ![]()
5: ![]()
6: ![]()
calculate ![]()
7: Calculate ![]()
8: ![]()
9: ![]()
10: for (i = 0 to![]()
) ![]()
11: Solve the following Pr by the dual simplex method to obtain ![]()
12: end if
13: else Optimized ¬ true//![]()
is an optimal solution to P.
14: end if
15: end while
3. Problem Instances and CPLEX Preprocessing
Four sets of NNLPs used in [6] are considered to evaluate the performance of the developed algorithm. Each problem set contains five problem instances for 21 different density levels and for varying ratios of (m constraints)/(n variables) from 200 to 1. Each set contains 105 randomly generated NNLPs with various densities p ranging from 0.005 to 1. Randomly generated real numbers between 1 and 5, 1 and 10, 1 and 10 were assigned to the elements of A, b and c respectively. To avoid having a constraint in the form of an upper bound on a variable, each constraint is required to have at least two non-zero![]()
. For general LP, a problem set containing 105 randomly generated by Saito et al. [17] is compared with the dynamic approach of this paper. These LP problems contain 1000 variables (n) and 200,000 constraints (m), with various densities ranging from 0.005 to 1 and the randomly generated ![]()
ranging between −1 and −5 or between 1 and 5.
Two parameters that CPLEX uses for solving linear programming are PREIND (preprocessing pre-solve indicator) and PREDUAL (preprocessing dual). As described in [17] and [6] , when parameter setting PREIND = 1 (ON), the preprocessing pre-solver is enabled and both the number of variables and the number of constraints is reduced before any type of algorithm is used. By setting PREIND = 0 (OFF) the pre-solver routine in CPLEX is disabled. PREDUAL is the second preprocessing parameter in CPLEX. By setting parameter PREDUAL = 0 (ON) or −1 (OFF), CPLEX automatically selects whether to solve the dual of the original LP or not. Both are used with the default settings for the CPLEX primal simplex method, the CPLEX dual simplex method, and the CPLEX barrier method. Neither CPLEX pre-solver nor PREDUAL parameters were used in any part of the developed dynamic active-set methods for NNLPs and LPs.
4. Computational Experiments
The computations were performed on an Intel Core (TM) 2 Duo X9650 3.00 GHz with a Linux 64-bit operating system and 8 GB of RAM. The developed methods use IBM CPLEX 12.5 callable library to solve linear programming problems. The dynamic RAD and dynamic GRAD are compared with the previously developed COST RAD and COST GRAD, respectively, as well as VIOL, the CPLEX primal simplex method, the CPLEX dual simplex method, and the CPLEX barrier method.
4.1. Computational Results for NNLP
Table 1 illustrates the performance comparison between dynamic RAD method and the previously defined constraint selection technique COST RAD on Set 1 to Set 4 for various dimensions of the matrix A used in [6] . Both methods are compared with the CPLEX barrier method (interior point), the CPLEX primal simplex method, and the CPLEX dual simplex method. The worst performance occurs at m/n ratio of 200, where on average, dynamic
![]()
Table 1. Results from dynamic RAD and COST RAD for set 1-set 4, (random, NNLP aij = 1 − 5, bi = 1 − 10, cj = 1 − 10).
+Used CPLEX preprocessing parameters of presolve = off and predual = off. ++Average of 5 instances of LPs at each density.
RAD is 8% faster than COST RAD for densities less than 0.2 and 18% slower for densities above 0.2. When the density increases, dynamic RAD shows an increase in computation time more than that of COST RAD. On the other hand, for an m/n ratio of 20 the CPU times decrease with an increase in density. For higher densities above 0.01, dynamic RAD is more efficient and takes less computation times than COST RAD. On average, dynamic RAD is 10% more efficient than COST RAD. For an m/n ratio of 2 at densities higher than 0.009, the data show that COST RAD starts taking significantly more time than dynamic RAD. Dynamic RAD was 5.5% faster than COST RAD over all densities and 21% faster on average for densities above 0.5. For an m/n ratio of 1 with densities greater than 0.01, dynamic RAD is about 8% more efficient than COST RAD. On average, dynamic RAD is superior performance to COST RAD for problem sets 2, 3, and 4.
Table 2 from [6] is presented to provide an immediate comparison of the developed dynamic RAD method with the standard CPLEX solvers. A reporting limit of 3000 seconds was used. On average, the CPU times for dynamic RAD were faster than any of the CPLEX solvers across all densities and ratios. However, CPLEX barrier methods show smaller CPU times when ratio m/n = 20 and the density is less than or equal to 0.01.
![]()
Table 2. Results from the CPLEX primal, the dual simplex, and the barrier method for set 1-set 4, (random, NNLP aij = 1 − 5, bi = 1 − 10, cj = 1 − 10) [6] .
−−Used CPLEX preprocessing parameters of presolve = ON and predual = Auto; ++Average of 5 instances of LPs at each density; bRuns with CPU times > 3000 s are not reported.
4.2. Computational Results for LP
Table 3 shows computational results for the CPLEX primal simplex method, the dual simplex method, and the interior point barrier method for the general LP problem set used in [17] . CPU times for COST GRAD and VIOL using both the multi-cut technique and dynamic approaches are presented for comparison. Dynamic GRAD is stable over the range of densities. In addition, its performance is superior to multi-cut GRAD for every problem instance. Average CPU times for GRAD using multi-cut method and dynamic approach are 43.87 and 24.57 seconds, respectively, a 42% improvement in computation time. Average computation times for GRAD and VIOL using dynamic approach are 24.57 seconds vs. 33.82 seconds, respectively.
It should be noted that GRAD captures more information than VIOL in higher densities to discriminate between constraints. Interestingly, when the dynamic active-set is used for both GRAD and VIOL, their CPU
![]()
Table 3. Comparison of computation times of CPLEX solvers, GRAD, and VIOL using both dynamic active-set and multi- cut method on general LP problem set (random LP with 1000 variables and 200,000 constraints [17] ).
+Used CPLEX preprocessing parameters of presolve = off and predual = off. 1Tx ≤ M = 1010 was used as the bounding constraint; ++Average of 5 instances of LPs at each density; −−Used CPLEX preprocessing parameters of presolve = ON and predual = Auto.
times are significantly faster than the same metrics with the multi-cut method. GRAD using the multi-cut technique takes the longest computation time in comparison to others at higher densities. Unlike the proposed dynamic approach, the LP algorithm
COST
GRAD requires checking the signs of the nonzero ![]()
and therefore more computation time for higher densities. The efficiency of VIOL decreases significantly with increasing density. On average, dynamic GRAD is approximately 35 times faster than the CPLEX primal simplex, 21 times faster than the CPLEX dual simplex, and 17 times faster the CPLEX barrier linear programming solvers without preprocessing. The superior overall performance of GRAD using dynamic approach is apparent across all densities in general LP set.
For comparison purposes, Table 4 shows GRAD and VIOL computation times when a fixed number of violated constraints is added at each iteration. Adding a fixed number of constraints is examined for both GRAD and VIOL. At densities below 0.03, dynamic GRAD takes less CPU time than the fixed-cut approach. GRAD
![]()
Table 4. Comparison of computation times of GRAD using dynamic active-set and fixed cut method on general LP problem set (random LP with 1000 variables and 200,000 constraints [17] ).
+Used CPLEX preprocessing parameters of presolve = off and predual = off. 1Tx ≤ M = 1010 was used as the bounding constraint; ++Average of 5 instances of LPs at each density.
*Corresponding author.
5. Conclusion
In this paper, dynamic active-set methods have been proposed for both NNLPs and LPs. In particular, these new approaches were compared to existing methods for problems with various sizes and densities. On average, dynamic RAD shows superior performance over COST RAD for the NNLP problem sets 2, 3, and 4. In the LP problem set, dynamic GRAD significantly outperformed the COST GRAD as well as the CPLEX primal simplex and the dual simplex. In this LP problem set, however, the barrier solver did outperform all methods for densities up to 0.03. In addition, dynamic GRAD outperformed a dynamic version of VIOL, which was a standard method in column generation and decomposition methods.
NOTES
*Corresponding author.