An Improved Search Approach for Solving Non-Convex Mixed-Integer Non Linear Programming Problems

The nonlinear mathematical programming problem addressed in this paper has a structure characterized by a subset of variables restricted to assume discrete values, which are linear and separable from the continuous variables. The strategy of releasing nonbasic variables from their bounds, combined with the “active constraint” method, has been developed. This strategy is used to force the appropriate non-integer basic variables to move to their neighbourhood integer points. Successful implementation of these algorithms was achieved on various test problems.


Introduction
There are a lot of optimization problem involved integer or discrete variables which occur naturally in a nonlinear function. Mixed Integer Nonlinear Programming (MINLP) refers to these kind of optimization problems with continuous and discrete variables and nonlinearities in the objective function and constraints. This paper concerns with a special class of MINLP problem in which the discrete variables are linear and separable from the continuous variables. This problem is defined by the following model. There are many diverse and important applications of MINLPs. A small subset of these applications includes portfolio optimization [1], [2], block layout design in the manufacturing and service sectors [3], network design with queuing delay constraints [4], integrated design and control of chemical processes [5], drinking water distribution systems security [6]. Energy system [7], [8], and multiperiod supply chain problems [9]. The needs in such diverse areas have motivated research and development in MINLP solver technology, particularly in algorithms for handling large-scale, highly combinatorial and highly nonlinear problems.
Methods for solving MINLPs include innovative approaches and related techniques taken and extended from MIP, such as, Outer Approximation (OA) methods [10,11], Branch-and-Bound (B&B) [12,13,24], Integrating Successive Quadratic Programming (SQP) and B&B [15], Generalized Bender's Decomposition (GBD) [16], and Branch and Cut [17] for solving MINLPs have been discussed in the literature since the early 1980's. These approaches generally rely on the successive solutions of closely related NLP problems. For example, B&B starts out forming a pure continuous NLP problem by dropping the integrality requirements of the discrete variables (often called the relaxed MINLP or RMINLP). Moreover, each node of the emerging B&B tree represents a solution of the RMINLP with adjusted bounds on the discrete variables. A survey regarding to solving a nonconvex MINLP can be found in [18].
Whenever the MINLP problem become too large to be solved by exact method, most of researchers move to heuristic method. Heuristic approaches to solving MINLPs include Variable Neighbourhood Search [19], automatically tuned variable fixing strategies [20], feasible neighbourhood search [21], Feasibility Pump [23,27], heuristics based on Iterative Rounding [25]. [19] propose a MINLP heuristic called the Relaxed-Exact-Continuous-Integer Problem Exploration (RECIPE) algorithm. The algorithm puts together a global search phase based on Variable Neighbourhood Search [19] and a local search phase based on a MINLP heuristic. Masmoudi et al. (017) propose a hybrid of genetic algorithm for solving MINLP. In heuristic approaches, however, one of the main algorithmic difficulties connected to MINLPs is to find a feasible solution.
Due to the fact that the functions in MINLPs are not smooth, therefore in this paper we use a direct search method, known as unconstrained optimization techniques that do not explicitly use derivatives. More information regarding to direct search method in optimization can be found in (Vayasis, 1991).
In this paper we address a strategy of releasing nonbasic variables from their bounds, combined with the "active constrained" method and the notion of superbasic for efficiently tackling a particular class of MINLP problems.
The rest of this paper is organized as follows. In Section 2 we give a brief notion of neighbourhood search. The basic approach of the proposed method is presented in Section 3. How to derive the proposed method is given in Section 4. The algorithm is presented in Section 5. Section 6 addresses a computational experience. The conclusions can be found in Section 7.

Neighbourhood Search
It should be noted that, generally, in integer programming the reduced gradient vector, which is normally used to detect an optimality condition, is not available, even though the problems are convex. Thus we need to impose a certain condition for the local testing search procedure in order to assure that we have obtained the "best" suboptimal integer feasible solution.
Scarf [18] has proposed a quantity test to replace the pricing test for optimality in the integer programming problem. The test is conducted by a search through the neighbors of a proposed feasible point to see whether a nearby point is also feasible and yields an improvement to the objective function.
Let Let [ ] be the integer feasible point which satisfies the above conditions. We could then say if [ ] + 1 ∈ ([ ]) implies that the point [ ] + 1 is either infeasible or yields an inferior value to the objective function obtained with respect to [ ]. In this case [ ] is said to be an "optimal" integer feasible solution to the integer programming problem. Obviously, in our case, a neigbourhood search is conducted through proposed feasible points such that the integer feasible solution would be at the least distance from the optimal continuous solution.

Solving Nonlinear Programming Problem
It has been our experience that many nonlinear programming problems are inordinately large because they are attempting to approximate, by piecewise linearization, what is essentially a nonlinear problem. It also appears that many real-life problems are such that only a small percentage of the variables are involved nonlinearly in the objective function. Therefore we should be able to write a class of nonlinear programming problems with the following form: Due to the separable property, w could partition x into a linear portion x L and a nonlinear portion x N . The partition can be expressed as where x N are called the nonlinear variables. Note that A and c operate on all variables x. In some cases the part of c T x involving x N may be incorporated into f(x N ); in other cases c may be zero. We assume that the function f(x N ) is continuously differentiable in the feasible region, with gradient  f(x N ) = g(x N ), and we assume that both f and g can be computed at any feasible point x N .
In essence the algorithm is an extension of the revised simplex method . To use some of the associated terminology, it might be described as an extension which permits more than m variables to be basic. Partitioning x and F(x) into linear and nonlinear terms is of considerable practical importance; for descriptive purposes, however, it is convenient to denote F(x) and F(x) simply by f(x) and g(x).
With a few conventional exceptions, we use upper-case letters for matrices, lower-case for vectors and Greek lower-case for scalars. The quantity  > 0 represents the precision of floating-point arithmetic.

Basic of The Method
Basic solutions; justification for standard form Before proceeding to the nonlinear problem, we need to introduce some linear programming background. In particular, equations (1)  (3) with f(x N ) = 0 are the standard form for stating linear programs and for solving them in practical implementations of the simplex method. A basic solution is characterized by having at most m "basic" variables lying between their bounds while the remaining n  m "nonbasic" variables are equal to one bound or other. An associated square basis matrix B is drawn from the columns of the constraint matrix A, and as the simplex method proceeds the columns of B are replaced one at a time.
Note that it is assumed in the standard form that A and x contain a full identity matrix and a full set of slack variables, respectively. (General equality and inequality constraints are accommodated by placing suitable upper and lower bounds on the slacks.) There are many practical reasons for retaining the standard form here. Full justification would require much implementation background, but very briefly, a sparse triangular factorization of B can be maintained more easily if columns (but not rows) of B are altered. Further, although the total number of variables is apparently greater, it is very easy to take advantage of the unit vectors associated with slacks, whenever B is re-factorized.

A. Superbasic variables
One virtue of the concept of basic solutions is the emphasis thereby given to the upper and lower bounds, l  x  u. It is misleading to regard these as "sparse constraints"; more importantly, they serve directly to eliminate a large percentage of the variables. The simplex method is therefore free to focus its attention on transforming (factorizing) just B, rather than the whole of A.
With nonlinear problems we cannot expect an optimal point to be a basic solution. However, if the number of nonlinear variables is small, it seems reasonable to suppose that an optimal solution will be "nearly" basic. Thus, as a simple generalization we introduce the notion of superbasic variables and partition the set of general constraints (2) as follows: (7) The matrix B is square and nonsingular as in the simplex method, S is m  s with 0  s  nm, and N is the remaining columns of A. The associated variables xB, xs, xN are called the basics, superbasics and nonbasics respectively. Both basics and superbasics are free to vary between their bounds. The name is chosen to highlight the superbasics' role as "driving force"; they may be moved in any direction at all (preferably one that improves the objective value), and the basics are then obliged to change in a definite way to maintain feasibility with respect to the constraints Ax = b.
Our expectation that solutions will often be "nearly" basic is confirmed by the following theorem: Theorem 1. Suppose a nonlinear program has t variables occurring nonlinearly (in either the objective or the constraints). An optimal solution exists in which the number of superbasic variables s satisfies s  t.

Derivation of method
We will assume that f(x) can be expanded in a Taylor's series with remainder of second order: where 0    1, and G(x + x) is the Hessian matrix of second partial derivatives evaluated at some point between x and x + x. Note that G is a constant matrix if f(x) is a quadratic function. Let us partition x and g(x) corresponding to the partitioning of A. If f(x) were truly quadratic, we could obtain a constrained stationary point at x + x by requiring two properties of the step x: i.e., the step remains on the surface given by the intersection of the active constraints.
i.e., the gradient at x + x (given by the left hand-side of (10)) is orthogonal to the surface of active constraints and is therefore expressible as a linear combination of the active constraint normals. For a more general function f(x), the step x may not lead directly to a stationary point, but we shall use Properties 1 and 2 to determine a feasible descent direction.
From (9) we have: and The following matrix can be used to simplify Eq. (10) 00 0 0 0 First it provides an expression for estimates of the Lagrange multipliers for the general constraints: (15) Note that when 0 s  x (which will mean x is stationary) we have As we observed that  is analogous to the pricing vector  in the revised simplex method. (From now on we shall denote the solution of (16) by .) Next we have from (10) x h (19) where 0 The form of equation (16)

Summary
We considered a class of algorithms in which the search direction along the surface of active constraints is characterized as being in the range of a matrix Z which is orthogonal to the matrix of constraint normals. Thus, if Apart from having full column rank, eq. (22) is (algebraically) the only constraint on Z and thus Z may take several forms. The particular Z corresponding to our own procedure is of the form This is a convenient representation which we will refer to for exposition purposes in later sections, but we emphasize that computationally we work only with S and a triangular (LU) factorization of B. The matrix Z itself is never computed.

The Algorithm
After solving the relaxed problem, the procedure for searching a suboptimal but integer-feasible solution from an optimal continuous solution can be described as follows. Let Otherwise go to next non-integer nonbasic or superbasic (if available). Eventually the column * is to be increased form LB or decreased from UB. If none go to next * . Step 4. Calculate * = −1 * i.e. solve * = * for * Step 5. Ratio test; there would be three possibilities for the basic variables in order to stay feasible due to the releasing of nonbasic * from its bounds. , describes a decision problem faced by a manufacturer of an electrical wire in order to upgrade its machinery. At present some of the machines are nearly worn out or obsolete, so that running costs are becoming unacceptably high. A shift in technology has caused an expansion in demand for the type of wire this manufacturer produces. However, the present plant cannot meet the demand projected for the near future.
There are two types of wire produced by this manufacturer: bare wire and covered wire. Each of these two products has two different sizes. In terms of the machines availability, the manufacturer has various options, including retaining and modifying some existing machines and purchasing new machines. Therefore it is desirable to draw a recommended configuration of the plant to meet the projected demand as well as to schedule which machine should produce each type of product such that the total running cost is minimized. The problem can be formulated mathematically as follows. The objective of the problem is to minimize the total running cost (in $1000's) 1  (COVRD2) Next is to make certain that the number of machines available for each option provides a production capacity that covers the production needs. Hence, the number of machine of type i is required for the production of all wire sizes.

Solution to the Problem
The plant upgrading problem mentioned in the previous section is a MILP problem. It has 11 constraints, 10 continuous variables and 5 integer variables. The relaxed problem is solved on PC optimization software. We then integerized the continuous result using our proposed algorithm. The continuous results as well as the integer results of the plant upgrading problem are presented in Table 1. The total computational time for both continuous and integer results is 2.5 seconds.  strategy of releasing nonbasic variable from their bounds, combined with the "active constraint" method and the notion of superbasic. After solving a problem by ignoring the integrality requirements, this strategy is used to force the appropriate non-integer basic variables to move to their neighborhoods integer points.
A study of the criteria for choosing a nonbasic variable to work with in the integerizing strategy has also been made. The number of integerizing steps would be finite if the number of integer variables contained in the problem is finite. However, it should be noted that the computational time for the integerizing process does not necessarily depend on the number of integer variables, since many of the integer variables may have an integer value at the continuous optimal solution.
The new direct search method has been shown to be successful on a range of problems, while not always able to achieve global optimality. In a number of cases to obtain the suboptimal point is acceptable, since the exponential complexity of the combinatorial problems in general precludes branchand-bound, except on small to medium problems.
Computational testing of the procedure presented this paper has demonstrated that it is a promised approach for large problems.