- Email: [email protected]
A soft approach for hard continuous optimization
European Journal of Operational Research 173 (2006) 18–29 www.elsevier.com/locate/ejor
Continuous Optimization
A soft approach for hard continuous optimization Chunhui Xu a
a,*
, Peggy Ng
b
Department of Management Information Science, Chiba Institute of Technology, Chiba 275-0016, Japan b School of Administrative Studies, York University, 4700 Keele Street, Toronto, Canada M3J 1P3 Received 12 November 2003; accepted 4 January 2005 Available online 25 February 2005
Abstract This paper is to introduce a soft approach for solving continuous optimizations models where seeking an optimal solution is theoretically or practically impossible. We first review methods for solving continuous optimization models, and argue that only a few optimization models with some good structure are solved. To solve a larger class of optimization problems, we suggest a soft approach by softening the goal in solving a model, and propose a two-stage process for implementing the soft approach. Furthermore, we offer an algorithm for solving optimization models with a convex feasible set, and verify the validity of the soft approach with numerical experiments. 2005 Elsevier B.V. All rights reserved. Keywords: Optimization; Optimal solution; Soft approach; Sampling
1. Introduction Optimization is on how to choose an alternative so as to optimize certain objective function. Various optimization problems can be simply formulated by the following model: min
*
J ðxÞ : x 2 X f ;
Corresponding author. Tel.: +81 47 478 0377; fax: +81 47 478 0354. E-mail address: [email protected] (C. Xu).
where x 2 Rn is a vector of decision variables, J(x) is an objective function, and X f is a set of feasible solutions, which is a huge or infinite set. Solving an optimization model traditionally means seeking an optimal solution x* in X f. Beginning from 1940s, much research had been done on solving models this way, and many successful stories have been reported, such as the simplex method for linear programming, hill-climbing algorithms for nonlinear programming. Section 2 gives a brief overview on these methods, and argues that only models with some good analytical structures have been solved.
0377-2217/$ - see front matter 2005 Elsevier B.V. All rights reserved. doi:10.1016/j.ejor.2005.01.004
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
When an optimization model is not of the good structures, how to solve the model has been a challenging problem. Section 3 introduces a new approach for optimization analysis, named a soft approach. After making clear the ideas in the soft approach, a two-stage process for implementing the soft approach is suggested. Section 4 offers an algorithm for solving continuous optimization models where the feasible set is a convex set. To illustrate the proposed soft approach, Section 5 applies the soft approach to solving two optimization problems. Section 6 concludes this paper and presents issues for further study.
2. Overview of methods for solving continuous models Let the feasible set in continuous optimization problems be ( ) n Y f X ¼ x2 X i jgi ðxÞ 6 0; i ¼ 1; . . . ; m ; ð1Þ i¼1
where Xi = [ai, bi] is the domain of xi, gi (x) is a constraint function (i = 1,. . ., m). The traditional linear programming, nonlinear programming, and the recent global optimization are the three branches on solving continuous optimization problems. Linear programming is an optimization model where the objective function and all constraint functions in the model are linear. The simplex method proposed by Danzig may be the most successful story in optimization, see the classic book by Danzig (1963) for more about LP. No more explanation on LP is necessary here, it is believed that LP is a solved model. Nonlinear programming is an optimization model where one or more nonlinear functions are included in the model. Many search methods, such as the NewtonÕs method and the conjugate gradient method, had been proposed for seeking an optimal solution. You can find a lot of examples where these algorithms work well, but their limitations are also obvious. First, they require all functions in the model to be differentiable. Second, they require the objective
19
function to be convex, otherwise, they only lead to stationary points or locally optimal solutions. Seeking an optimal solution of a nonconvex or/and nondifferentiable model has been a difficult problem. See Lasdon (1970), Luenberger (1984) and Fletcher (1987) for details about nonlinear programming. Global optimization is a relatively new field which aims at seeking an optimal solution in nonlinear optimization models. The main stream in global optimization exploits analytical properties of the model to generate a deterministic sequence of points converging to a global optimal solution. Two analytical properties, convexity and monotonicity, have been successfully exploited, giving rise to two research trends in global optimization: d.c. optimization and monotonic optimization. The d.c. optimization deals with models described by differences of convex functions, while monotonic optimization deals with models described by means of functions monotonically increasing or decreasing along rays, hence the proposed algorithms still depend on these analytical properties of models. Refer to Horst et al. (1995), Horst and Pardalos (1995), for more about global optimization. Therefore, we can see that only models with some good structures have been solved. For the general nonlinear model, how to find out an optimal solution remains to be a challenging problem. Remark 1. Other methods for solving optimization problems have been proposed that do not necessarily aim at seeking an optimal solution with certainty, such as the stochastic approaches and the metaheuristics. Stochastic approaches generally require little or no additional assumptions on the structure of optimization models, but the outcome they lead to is random, although some provide a probabilistic convergence guarantee, refer to Pardalos et al. (2000) for a review of three kinds of stochastic methods. These approaches are not practically usable since they can not quantify what a result they will produce within a given time period. Metaheuristics are methods that are often based on processes observed in physics or biology. Three main categories of them are: simulated annealing,
20
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
tabu search, and genetic algorithms. We do not dwell on metaheuristics here, but refer readers to Laarhoven and Aarts (1987) for simulated annealing, Glover and Lagunna (1997) for Tabu search, and Goldberg (1989) for genetic algorithms. Many experiments showed that metaheuristics were very successful in solving hard combinatorial optimization problems, but they are heuristic in nature, and they can not quantify the result they lead to. The soft approach to be introduced in this paper differs essentially from the stochastic approaches and the metaheuristics because we can quantify the result obtained, which will be explained in the following sections.
3. A soft approach As stated in the previous section, only a small part of optimization models has been solved. To make matters worse, some ‘‘solved’’ models may not be really solved if taking some practical considerations into account. For instance, when there is a time limit for computation, algorithms which can not lead to an optimal solution within the given time limit will be helpless, because an algorithm converging to an optimal solution in 10 minutes may produce a very bad result when terminated in 8 minutes. Regretfully, most traditional algorithms are not capable to handle time limit requirements. For solving a larger class of optimization models, and also for handling the time limit requirement, this section introduces a soft approach, which does not require any assumption about the objective function, and can handle requirements on computation time. We first explain the ideas in the soft resolution approach, and then suggest a two-stage process for implementing it.
optimization, which was proposed for solving complicated stochastic optimization models, refer to Lau and Ho (1997), and Lee et al. (1999) for details on ordinal optimization. This goal softening method softens the goal of solving an optimization model in two dimensions: 1. Instead of looking for an optimal solution, it seeks a good enough solution. To define a good enough solution, we donÕt use the cardinal performance but the order of a solution in the feasible set. For instance, the top k% solutions may be taken as good enough solutions, where k is set by the requirements in solving a model. 2. Instead of requiring the solution to be an optimal solution for sure, it accepts a solution if it is high likely to be a good enough solution. When the probability with which an alternative is a good enough solution is larger than some value, say 95%, this alternative can be taken as a solution. In other words, the soft approach agrees that a model is ‘‘solved’’ if an alternative which is high likely a good enough solution is obtained. The notions of good enough and high probability are set by the requirements in solving a model. To focus on the main issues, letÕs take the top 1% solutions as good enough, and 95% as high probability. Therefore, when a solution is founded to be in the top 1% with a probability bigger than 95%, the model is solved. The good enough solutions usually form a subset of X f, seeking any point of a set is much easier than seeking exactly one point, it is easy to understand that ‘‘solving’’ an optimization model under the soft approach should be easier. In addition, the soft approach allows failure, i.e., we do not require the solution found to be a good enough solution for sure, it is reasonable to believe that a large class of models are expected to be solvable under the soft approach.
3.1. Ideas in the soft approach 3.2. A two-stage process The soft approach softens the goal in ‘‘solving’’ an optimization model. Among several goal softening methods, we use the one suggested in ordinal
We suggest to implement the soft approach in the following two stages.
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
• Stage 1: Sample the feasible set X f to generate a finite subset S. As for the sample set S, we require the probability with which S contains at least one good enough solution to be high. Let G X f be a set of good enough solutions, say the top 1% solutions, we require the following condition to be satisfied: pfjS \ Gj P 1g P 95%:
ð2Þ
Consequently, the best alternative in S will be a good enough solution with a probability not smaller than 95%. • Stage 2: Select the best alternative x* from the sample set S. Then x* is the solution we are looking for. This two-stage process is illustrated in Fig. 1. Since S is a finite set, there is no conceptual difficulty in Stage 2. However, Stage 2 will be time consuming if S is too big. Fortunately, we will show below that a small sample set is enough to satisfy condition (2), no need to generate a big sample set, so Stage 2 will not take much computation time. Suppose that we generate S by uniformly taking 1000 samples from X f, then the probability that none of the good enough solutions is included in S is p = (1 1%)1000 < 0.1%, consequently, p{jS \ Gj P 1} = 1 p > 99.9%. That is, a sample set with 1000 uniform samples can satisfy (2). In other words, if we generate 1000 uniform samples from X f in Stage 1, and choose the best sample in Stage 2, then the best sample is one of the top 1% solutions with a probability not smaller than 99%, so the model is solved.
21
3.3. Features of the soft approach From the two-stage process, we can easily see that the soft approach has the following two distinguished features. • No assumption about the objective function is needed. The objective function is used only in Stage 2. Since S is a finite and small set, say with 1000 samples, we can evaluate all samples in S to get the best sample, so no assumption about the objective function is needed. • Time for solving a model is controllable. We can estimate the computation time needed for solving an optimization model given the requirements in solving the model, and quantify the quality of the solution obtained within a given time limit. Let t1 be the time for taking one feasible sample, and t2 be the time for evaluating one sample, then the time for solving a model, denoted by T, is equal to T ¼ jSjðt1 þ t2 Þ;
ð3Þ
where jSj is the number of samples taken in Stage 1, which should be determined by the requirements in solving an optimization model and the sampling method used in Stage 1. Suppose it is required that the solution should be one of the top k% alternatives with a probability not smaller than q%, and uniform sampling is used in Stage 1, then jSj can be obtained by solving the following equation: jSj
O O
O O O O
O O O O O
S ** * * * * * ** O * G O
O
G: Set of good enough solutions S: Set of samples O Sampling
O
O Selecting O O O
*
** * S o: Feasible solution *: Good enough solution
Fig. 1. A two-stage process.
pfjS \ Gj P 1g ¼ 1 ð1 k%Þ
¼ q%:
ð4Þ
T is known by substituting jSj into (3). So we can estimate the time needed for solving a model when requirements of solving a model are given. For instance, suppose the solution obtained should be one top 1% alternative with a probability bigger than 99.99%, then we know from (4) that the minimal number of uniform samples we need to take in Stage 1 is logð199:99%Þ ¼ 917. If generatlogð11%Þ ing one uniform sample takes 103 seconds and
22
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
calculating the value of objective function at each sample point takes 104 seconds with a computer, then total computation time needed to solve the problem with this computer is approximately 917 * (103 + 104) 1 second. On the other hand, we can quantify the quality of the solution obtained given a time limit for solving a model. Let T0 be the time limit given, then the number of samples we can take in Stage 1 is jSj ¼
T0 : t1 þ t2
ð5Þ
When Stage 1 takes samples with some uniform sampling method, substituting (5) into the following equation: jSj
pfjS \ Gj P 1g ¼ 1 ð1 k%Þ ;
ð6Þ
we can know the probability with which the solution is one of the top k% solutions. For instance, suppose we are given 2 seconds for solving a model, and it takes a computer 103 seconds for generating one uniform sample and 104 seconds for calculating the value of objective function at each sample point, then the maximum number of samples we can take in Stage 2 1 is 103 þ10 4 ¼ 1818. If the good enough solutions are defined as the top 0.1% points, then the soft approach can generate a good enough solution in 2 seconds with a probability of 1(10.1%)1818 = 83.78%, thus the quality of the solution obtained within 2 seconds is quantified. Remark 2. Time needed for taking one feasible sample depends on three factors: the sampling method in Stage 1, the feasible set of a model and computer speed. Given these three factors, we can know t1 by doing a simple sampling test. Time needed for evaluating one sample depends on the speed of a computer and the objective function of a model, we anticipate that t2 will be much smaller than t1, so time needed in solving an optimization model mainly depends on t1. Computing time controllability may be most valuable feature of the soft approach, cause other optimization methods, including traditional optimization methods, stochastic proaches and metaheuristics, can not handle
the bethe apthe
requirement on computing time, and there are many situations where models need to be solved within some time limit. The soft approach not only gives a solution within a give time limit, but also tells the quality of the solution obtained.
4. An algorithm for solving continuous models This section presents an algorithm for solving continuous models following the two-stage process. To complete Stage 1, we need to decide the number of necessary samples and a method for taking samples from X f. Suppose that the requirement for solving a model is: seeking one of the top k% solutions with a probability not small then q%. Then the sample set S needs to satisfy the following condition: pfjS \ Gj P 1g P q%:
ð7Þ
Suppose that no information is available about the distribution of the top k% solutions, and Stage 1 uses some uniform sampling method. Then the smallest number of necessary samples is determined by solving the following equation with jSj as the variable: 1 ð1 k%ÞjSj ¼ q%:
ð8Þ
We next consider methods for uniform sampling. Let ð9Þ X f ¼ fx 2 X jgj ðxÞ 6 0; j ¼ 1; . . . ; mg; Qn where X ¼ i¼1 X i is the domain of x. Without loosing generality, we assume that points satisfying all constraints are also in the domain.1 There proposed two approaches for uniformly sampling a set described by (9): the rejection approach or the random generating approach. The rejection approach takes samples from X uniformly and accepts samples if they are also in X f and reject otherwise. This method applies to uniformly sample any bounded sets. However, when the feasible set X f is just a very small part of the domain, a lot of points will be rejected be1
This means that X1 = {x 2 Rnjgi(x) 6 0, i = 1,. . ., m} X holds. When X1XX, we can add constraints like ai 6 xi 6 bi to X1 so that X1 X will hold.
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
fore a feasible point is founded, Stage 1 will be time consuming due to the low sampling efficiency. The random generating approach suggests generating uniform samples within the feasible set, its efficiency is high if a proper generating method can be worked out. Smith (1984) proposed a general mixing method for generating uniform samples in a bounded set, which is stated below. • 1. 2. 3.
A general mixing method Find a start point x0 2 X f, let i = 0. Generate a random direction d in Rn. Form a line set L = X f \ {xjx = xi + kd}, where k is a real scalar. 4. Generate a random point xi+1 uniformly distributed over L. 5. Set i to i + 1, go back to Step 2 till necessary samples are obtained.
This general mixing method is illustrated in Fig. 2. This mixing method is good in generality, but it did not tell how to determine k in Step 3. Berbee et al. (1987) proposed a method called ‘‘hit-and-run’’ algorithm for uniformly sampling a simplex. The main difference with the above algorithm is in Step 3, the hit-and-run algorithm forms the line L by connecting two hit points in the boundary of a simplex, which are decided by hitting the boundary of a simplex along directions d and d. We employ this idea to develop an algorithm for uniformly sampling a convex set. The main idea of our algorithm is as follows: Since X f is convex, hitting the boundary of X f from an interior point along one direction will produce only one feasible hit point, and the line seg-
ment connecting any two feasible hit points is still in X f. For a uniformly generated direction d, we can get two hit points in boundary, one in direction d and the other in direction d, then the point uniformly distributed on the line segment connecting the two hit points is a uniform sample. The algorithm is stated in the following. • A uniform sampling algorithm for convex sets 1. Find an interior point x0 2 X f uniformly as a start point, let i = 0. x0 can be founded by using the rejection method while uniformly sampling X. 2. Select a direction d uniformly over a unit sphere D = {d 2 Rnjjjdjj = 1}.2 3. Find the two hit points by hitting the boundary of X f from standing point xi along directions d and d. Denote Pi, Qi: hit points in the boundary of X f along directions d and d, respectively, k1: distance between xi and Pi, k2: distance between xi and Qi. Then we have P i ¼ xi þ k1 d;
d xi
x
i+1
Qi ¼ xi k2 d: 1
ð10Þ
2
We can calculate k and k as follows: Find hit points in each constraint surface gj (x) = 0, wthe distance between xi and the nearest hit point in the constraint surface along direction d is k1, while the distance along direction d is k2. For each j 2 {1,. . .,m}, denote k1j : k2j :
L
23
distance between xi and the nearest hit point in constraint surface gj (x) = 0 along direction d, distance between xi and the nearest hit point in constraint surface gj (x) = 0 along direction d.
Then we know that k1j is the smallest positive solution of the following equation:
Xf Fig. 2. A general mixing process.
2 Taking uniform samples from nice regions like a hyper cube or a hypersphere is easy, refer to Rubinstein (1982).
24
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
gj ðxi þ kdÞ ¼ 0;
ð11Þ
while k2j is the smallest positive solution of gj ðxi kdÞ ¼ 0;
ð12Þ
When (11) or/and (12) do not have a position solution, set k1j or/and k2j to 1, thus we have k1 ¼ minfk11 ; ; k1m g; k2 ¼ minfk21 ; ; k2m g:
For the case where all constraints are linear, finding the hit points is quite easy. Refer to Berbee et al. (1987) for details. When necessary samples are obtained using the above sampling algorithm, Stage 2 selects the best alternative by evaluating all samples in S, leading to one alternative which meets the requirements in solving the model.
ð13Þ
i+1
uniformly over the line segment 4. Choose x connecting the two hit points: xiþ1 ¼ xi þ lðP i Qi Þ;
ð14Þ
where l is a random variable distributed uniformly over [0,1]. 5. Set i to i + 1, go back to Step 2 till necessary samples are obtained. This algorithm is illustrated in Fig. 3. Remark 3. To ensure getting two hit points in each iteration, it is necessary to hit the boundary from an interior point, two strategies for doing this are: (1) Always hit the boundary from the start point x0, (2) If a hit point in the boundary is chosen in Step 4 (l = 0 or 1), do not change standing point in next hitting. When there is not an interior point in X f, refer to Xu and Ho (2002) for a uniform sampling method.
5. Numerical examples We will use the soft approach to solve two optimization problems for illustrating the validity of the soft approach. Example 1. Find the minimum of the following function in interval [2, 2]: ( xj sin ð1xÞj ; x 6¼ 0; x2 þ1 f ðxÞ ¼ ð15Þ 0; x ¼ 0: The figure of function (15) is shown in Fig. 4. This is a sort of hard optimization problem since the objective function is not differentiable at many points, and there are many local minimums. For testing the soft approach, we solved this problem with some other method, the minimum is reached at x0 = 0.7374 and f (x0) = 0.4667. To see how well the soft approach solves this model, suppose that the set of good enough
Remark 4. Seeking hit points is the main computation load in this sampling algorithm, which requires to solve Eqs. (11) and (12) for each constraint function with k as the variable. Since these equations contain only one variable k, there should be no essential difficulty in solving these equations.
0.5 0.4 0.3 0.2
f(x)
0.1 0 -0.1
P
1
d
Q
i
i
x i+1
2
d
-0.2 -0.3 -0.4
xi X
f
Fig. 3. A uniform sampling algorithm for convex sets.
-0.5 -2
-1.5
-1
-0.5
0
0.5
1
1.5
x
Fig. 4. The objective function in Example 1.
2
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
G fx 2 Rjjx x0 j 6 0:04g:
ð16Þ
f
Since X is an interval [2, 2], uniformly sampling X f is very easy, we do not need to use the sampling algorithm proposed in Section 4. After getting 1000 uniform samples in Stage 1, we choose the best sample as the solution. We solved this optimization problem 50 times with the soft approach, and got one good enough point each time. Thus we have illustrated the effectiveness of the soft approach in this example. Table 1 below summarizes the results produced, wherein x* is the solution obtained, and f (x*) is the value of the objective function at point x*. For showing the results clearly, we plot in Fig. 5 the objective values of the best samples in the 50 experiments.
-0.4666 -0.4666 -0.4666 -0.4667
f(x*)
solutions G is defined as the top 1% points in X f = [2, 2], and the solution produced should be a good enough solution with a probability bigger than 99.99%. When taking 1000 uniform samples from X f, the possibility that at least one point in G is included in the samples is p(jG \ X fj P 1) = 10.991000 = 99.996%, so taking 1000 uniform samples in Stage 1 is enough. To check if the solution obtained from using the soft approach is in G, we need to determine G and do some computation experiments. Let the length of an interval be its metric, then the metric of G, denoted by jGj, is 1% * jX fj = 0.04. Fig. 4 shows that the objective function is nearly symmetric in the neighborhood of x0, we can take the points in the neighbor of x0 as an approximation of G, that is,
25
-0.4667 -0.4667 -0.4667 -0.4667 -0.4667 0
5
10
15
20
25
30
35
40
45
Fig. 5. Objective values of the best samples in 50 experiments (Example 1).
We can also examine the degree of approximation in term of value since we know the optimal solution. Let the degree of approximation be measured by f ðx Þ f ðx0 Þ ; d ¼ f ðx0 Þ where x* is a solution obtained from using the soft approach, x0 is the optimal solution. We calculated d in each of the 50 experiments, the worst approximation turned out to be 0.0122%, so the soft approach leaded to a good enough solution even if in term of value. Example 2. Find the minimum of the following function: f ðx1 ; x2 Þ ¼
j sinð5px1 Þ þ sinð5px2 Þj e10ðx1 0:25Þ
2
þ10ðx2 0:5Þ2
ð17Þ
;
Table 1 Computation results from solving Example 1 x* f (x*) x* f (x*) x* f (x*) x* f (x*) x* f (x*)
0.7384 0.4667 0.7305 0.4666 0.7380 0.4667 0.7376 0.4667 0.7369 0.4667
0.7380 0.4667 0.7353 0.4667 0.7374 0.4667 0.7376 0.4667 0.7412 0.4667
0.7370 0.4667 0.7393 0.4667 0.7362 0.4667 0.7375 0.4667 0.7363 0.4667
0.7368 0.4667 0.7364 0.4667 0.7374 0.4667 0.7396 0.4667 0.7409 0.4667
0.7386 0.4667 0.7386 0.4667 0.7380 0.4667 0.7392 0.4667 0.7384 0.4667
50
Sequence of experiments
0.7378 0.4667 0.7384 0.4667 0.7385 0.4667 0.7407 0.4667 0.7361 0.4667
0.7334 0.4667 0.7363 0.4667 0.7384 0.4667 0.7339 0.4667 0.7368 0.4667
0.7324 0.4667 0.7341 0.4667 0.7368 0.4667 0.7302 0.4666 0.7404 0.4667
0.7444 0.4666 0.7380 0.4667 0.7382 0.4667 0.7379 0.4667 0.7387 0.4667
0.7378 0.4667 0.7394 0.4667 0.7402 0.4667 0.7380 0.4667 0.7354 0.4667
26
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
in the feasible set defined by 2
X f ¼ fðx1 ; x2 Þ 2 ½0; 1 j0 6 1 x1 x2 6 1; x1 þ 2x2 P 1g:
ð18Þ
The figure of the objective function (17) is shown in Fig. 6, and X f is shown in Fig. 7. This is a sort of hard problem since the object function is not convex. For testing the soft approach, we solved this problem with some other method, the minimum is reached at (0.293, 0.5) and f (x0) = 1.9574, that is, x0 = (0.293, 0.5), f (x0) = 1.9574. Similar as in Example 1, let the set of good enough solutions be the top 1% points in X f, we
0
f(x1,x2)
-0.5 -1 -1.5
are satisfied with obtaining one good enough solution with a probability bigger than 99.99%, hence we need only taking 1000 uniform samples from X f in Stage 1. Fig. 7 shows that X f accounts for only 1/4 of the domain [0,1]2, we use the sampling algorithm proposed in Section 4 for sampling X f efficiently. Fig. 8 shows the 1000 samples taken uniformly from X f with the proposed algorithm. After getting 1000 uniform samples, choose the best sample as the solution. In order to check if the soft approach produces a good enough solution with a high probability, we need to determine G and do computation experiments. Let the area of a set in R2 be its metric, then jX fj = 13/4 = 0.25, and jGj = jX fj * 1% = 0.0025. Fig. 6 shows that the objective function is nearly symmetric in the neighborhood of x0, we take the circle around x0 as an approximation of G, that is, G {x 2 R2jjjxx0jj 6 r}. Since jGj = 0.0025, r is determined from solving equation pr2 = 0.0025, that is, r = 0.00282, hence we have G fx 2 R2 jjjx x0 jj 6 0:0282g:
-2 1 0.8
1
0.6 0.4
x2
0.2 0
0.2
0
0.4
0.6
0.8
x1
ð19Þ
We solved this optimization problem 50 times with the soft approach, and got one good enough point each time. Thus we have illustrated the effectiveness of the soft approach in this example.
1
1
0.9
0.9
0.8
0.8
0.7
0.7
0.6
0.6
X2
X2
Fig. 6. The objective function in Example 2.
0.5
0.5
0.4
0.4
0.3
0.3
0.2
0.2
0.1
0.1
0
0
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
X1
Fig. 7. Feasible set in Example 2.
0.8
0.9
1
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
X1
Fig. 8. Uniform samples from X f.
0.8
0.9
1
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
27
Table 2 Computation results from solving Example 2 x1 x2 f (x*) x1 x2 f (x*) x1 x2 f (x*) x1 x2 f (x*) x1 x2 f (x*)
0.2886 0.4997 1.9546 0.2949 0.5008 1.9568 0.2882 0.4960 1.9519 0.3009 0.5031 1.9475 0.2898 0.4920 1.9471
0.2958 0.5084 1.9465 0.2861 0.4982 1.9501 0.2928 0.5062 1.9519 0.2978 0.5018 1.9537 0.2987 0.5086 1.9424
0.2867 0.4913 1.9413 0.2949 0.4909 1.9454 0.2979 0.4956 1.9514 0.2878 0.5045 1.9507 0.2917 0.4975 1.9563
0.2993 0.5045 1.9490 0.2897 0.4962 1.9539 0.2881 0.4927 1.9465 0.2866 0.4921 1.9429 0.2946 0.4965 1.9554
0.3013 0.5042 1.9455 0.2873 0.5016 1.9524 0.3029 0.4956 1.9411 0.2915 0.5061 1.9519 0.2874 0.4996 1.9530
Table 2 summarizes the results produced, wherein x* is the solution obtained, and f (x*) is the value of the objective function at point x*. For seeing the results clearly, we plot in Fig. 9 the objective values of the best samples in the 50 experiments. We can also examine the degree of approximation in term of value of since we know the optimal solution. We calculated the degree of approximation in each of the 50 experiments, the worst degree of approximation turned out to be 1.9752%, so the soft approach leaded to a good solution even if in term of value. -1.915 -1.92 -1.925 -1.93
f(x*)
-1.935 -1.94 -1.945 -1.95 -1.955 -1.96
0
5
10
15
20
25
30
35
40
45
50
Sequence of experiments
Fig. 9. Objective values of the best samples in 50 experiments (Example 2).
0.2933 0.5168 1.9182 0.2895 0.5013 1.9554 0.2885 0.4981 1.9541 0.2870 0.5003 1.9523 0.2898 0.4991 1.9559
0.3096 0.5005 1.9190 0.3015 0.5039 1.9451 0.2978 0.5059 1.9494 0.2946 0.5012 1.9569 0.2815 0.4986 1.9384
0.2854 0.5029 1.9481 0.2979 0.5001 1.9541 0.2907 0.5019 1.9561 0.3057 0.4954 1.9322 0.3062 0.5049 1.9299
0.2831 0.5008 1.9435 0.2954 0.5052 1.9529 0.2875 0.4920 1.9442 0.2815 0.4907 1.9263 0.2960 0.5008 1.9561
0.2973 0.5060 1.9498 0.2864 0.5022 1.9506 0.2982 0.5028 1.9526 0.2997 0.4874 1.9292 0.2933 0.4966 1.9558
Remark 5. We specified the good enough set in the two examples, this is only for checking if the result obtained is G. The soft approach does not require that we should specify G, and it is also not necessary to know where the good enough solutions are located in the feasible set in using the soft approach. We got the optimal solution with other optimization methods in the two examples, this is for specifying the good enough set. It is not necessary to know the optimal solution for quantifying the solution produced by the soft approach, because the soft approach does not use value in quantifying the solution produced. Computation time of the soft approach is determined by the time needed for generating necessary samples and calculating the values of the objective function at the sample points, which depends on the speed of computer used. We used a IBM ThinkPad computer (CPU: Pentium M Processor 1.6 GHz, Memory: 768MB) in solving the above model. It takes 1.5 * 104 seconds in generating one uniform sample with the proposed method, and 7 * 105 seconds in evaluating one sample, so solving this model takes 1000 * (1.5 * 104 + 7 * 105) = 0.22 seconds. To see what a solution the soft approach can produce within a time limit, we do some computing experiments by limiting the computing time to 2 seconds.
28
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
With the same computer, the maximum number samples we can take in Stage 1 is 2 ¼ 9091. Suppose the good enough 1:5104 þ7105 solutions are defined as the top 0.1% points, then the soft approach is expected to produce within 2 seconds one good enough solution with a probability of 1(10.1%)9091 = 99.99%. To verify this theoretical anticipation, we specify G first and then do 50 computing experiments. Since jGj = 0.1%jX fj = 0.00025 and G {x 2 R2jjjxx0jj 6 r}, we can know r from solving equation pr2 = 0.00025, that is, r = 0.0089, hence of
G fx 2 R2 jjjx x0 jj 6 0:0089g:
ð20Þ
We solve this optimization problem by taking 9091 uniform samples in Stage 1 and choosing the best sample as the solution in Stage 2. We repeated this solving process 50 times, and got one point in G defined by (20) each time, thus the soft approach did produce a good enough solution within a given time limit. The objective values of the best samples obtained from the 50 experiments are shown in Fig. 10. We omit the detailed experiment data for saving space. We calculated the degree of approximation in each of the 50 experiments, and the worst degree of approximation turned out to be 0.31%. Comparing with the worst degree of approximation
6. Conclusions This paper introduced a soft approach for solving complicated optimization models, and suggested a two-stage process for implementing the soft approach. An algorithm was given for solving optimization problems with a convex feasible set. The soft approach has two distinguished features that most traditional optimization methods do not have, i.e., no assumption about the objective function is needed, and the computation time is controllable. The soft approach differs essentially from stochastic and heuristic methods in the sense that the result obtained can be quantified. The soft approach does not intend to replace the optimum seeking methods, it does provide an alternative when optimum seeking fails. We expect that the soft approach will be helpful when no method is available for seeking an optimal solution, or when existing algorithms can not generate a good solution within a given time limit. Further research on the soft approach can be conducted in the following directions: (1) Efficient sampling methods. Sampling plays the central role in the soft approach, designing efficient sampling methods for different sets is the key for getting a better result. (2) Biased sampling methods. When there is some information about the distribution of good enough solutions, concentrating on the areas with a high density of good enough solutions in sampling will produce better results than uniform sampling, how to design such sampling methods is an interesting issue.
-1.951
-1.952
-1.953
-1.954
f(x*)
when taking 1000 samples, which was 1.9752% in our previous 50 experiments, we can see that taking more samples in Stage 1 improved the quality of solution even in term of value.
-1.955
-1.956
-1.957
-1.958
Acknowledgments 0
5
10
15
20
25
30
35
40
45
50
Sequence of experiments
Fig. 10. Objective values of the best samples in 50 experiments when taking 9091 samples.
Dr. C. Xu is grateful to Professor Y.C. Ho at Harvard University for his kind support and insightful guidance while doing research with him
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
at Harvard, and also for the wonderful ideas in his ordinal optimization, from which this study benefits a lot.
References Berbee, H.C.P., Boender, C.G.E., Kan, A.H.G.R., Scheffer, R.L., Smith, R.L., Telgen, J., 1987. Hit-and-run algorithms for the identification of nonredundant linear inequalities. Mathematical Programming 37 (2), 184–207. Danzig, G.B., 1963. Linear Programming and Extensions. Princeton University Press, Princeton, NJ. Fletcher, R., 1987. Practical Methods of Optimization. John Wiley & Sons, New York. Glover, F., Lagunna, M., 1997. Tabu Search. Kluwer Academic Publishers, Boston, MA. Goldberg, D., 1989. Genetic Algorithms in Search, Optimization, and Machine Learning. Addision-Wesley, San Mateo, CA. Horst, R., Pardalos, P.M., Thoai, N.V., 1995. Introduction to Global Optimization. Kluwer Academic Publishers, Dordrecht. Horst, R., Pardalos, P.M. (Eds.), 1995. Handbook of Global Optimization. Kluwer Academic Publishers, Dordrecht.
29
Laarhoven, P.J.M. van, Aarts, E.H., 1987. Simulated Annealing: Theory and Applications. Kluwer Academic Publishers, Dordrecht. Lasdon, L.S., 1970. Optimization Theory for Large Systems. Macmillan, London. Lau, T.W.E., Ho, Y.C., 1997. Universal alignment probabilities and subset selection for ordinal optimization. JOTA 93 (3), 455–489. Lee, L.H., Lau, T.W.E., Ho, Y.C., 1999. Explanation of goal softening in ordinal optimization. IEEE Transactions on AC 44 (1), 94–99. Luenberger, D.G., 1984. Introduction to Linear and Nonlinear Programming. Addison-Wesley, Reading, MA. Pardalos, P.M., Romeijn, H.E., Tuy, H., 2000. Recent developments and trends in global optimization. Journal of Computational and Applied Mathematics 124, 209–228. Rubinstein, R.Y., 1982. Generating random vectors uniformly distributed inside and on the surface of different regions. European Journal of Operational Research 10 (2), 205–209. Smith, R.L., 1984. Efficient Monte-Carlo procedures for generating points uniformly distributed over bounded regions. Operations Research 32 (6), 1296–1308. Xu, C., Ho, Y.C., 2002. Constrained ordinal optimization, Unpublished research report, Division of Engineering and Applied Sciences, Harvard University.
Continuous Optimization
A soft approach for hard continuous optimization Chunhui Xu a
a,*
, Peggy Ng
b
Department of Management Information Science, Chiba Institute of Technology, Chiba 275-0016, Japan b School of Administrative Studies, York University, 4700 Keele Street, Toronto, Canada M3J 1P3 Received 12 November 2003; accepted 4 January 2005 Available online 25 February 2005
Abstract This paper is to introduce a soft approach for solving continuous optimizations models where seeking an optimal solution is theoretically or practically impossible. We first review methods for solving continuous optimization models, and argue that only a few optimization models with some good structure are solved. To solve a larger class of optimization problems, we suggest a soft approach by softening the goal in solving a model, and propose a two-stage process for implementing the soft approach. Furthermore, we offer an algorithm for solving optimization models with a convex feasible set, and verify the validity of the soft approach with numerical experiments. 2005 Elsevier B.V. All rights reserved. Keywords: Optimization; Optimal solution; Soft approach; Sampling
1. Introduction Optimization is on how to choose an alternative so as to optimize certain objective function. Various optimization problems can be simply formulated by the following model: min
*
J ðxÞ : x 2 X f ;
Corresponding author. Tel.: +81 47 478 0377; fax: +81 47 478 0354. E-mail address: [email protected] (C. Xu).
where x 2 Rn is a vector of decision variables, J(x) is an objective function, and X f is a set of feasible solutions, which is a huge or infinite set. Solving an optimization model traditionally means seeking an optimal solution x* in X f. Beginning from 1940s, much research had been done on solving models this way, and many successful stories have been reported, such as the simplex method for linear programming, hill-climbing algorithms for nonlinear programming. Section 2 gives a brief overview on these methods, and argues that only models with some good analytical structures have been solved.
0377-2217/$ - see front matter 2005 Elsevier B.V. All rights reserved. doi:10.1016/j.ejor.2005.01.004
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
When an optimization model is not of the good structures, how to solve the model has been a challenging problem. Section 3 introduces a new approach for optimization analysis, named a soft approach. After making clear the ideas in the soft approach, a two-stage process for implementing the soft approach is suggested. Section 4 offers an algorithm for solving continuous optimization models where the feasible set is a convex set. To illustrate the proposed soft approach, Section 5 applies the soft approach to solving two optimization problems. Section 6 concludes this paper and presents issues for further study.
2. Overview of methods for solving continuous models Let the feasible set in continuous optimization problems be ( ) n Y f X ¼ x2 X i jgi ðxÞ 6 0; i ¼ 1; . . . ; m ; ð1Þ i¼1
where Xi = [ai, bi] is the domain of xi, gi (x) is a constraint function (i = 1,. . ., m). The traditional linear programming, nonlinear programming, and the recent global optimization are the three branches on solving continuous optimization problems. Linear programming is an optimization model where the objective function and all constraint functions in the model are linear. The simplex method proposed by Danzig may be the most successful story in optimization, see the classic book by Danzig (1963) for more about LP. No more explanation on LP is necessary here, it is believed that LP is a solved model. Nonlinear programming is an optimization model where one or more nonlinear functions are included in the model. Many search methods, such as the NewtonÕs method and the conjugate gradient method, had been proposed for seeking an optimal solution. You can find a lot of examples where these algorithms work well, but their limitations are also obvious. First, they require all functions in the model to be differentiable. Second, they require the objective
19
function to be convex, otherwise, they only lead to stationary points or locally optimal solutions. Seeking an optimal solution of a nonconvex or/and nondifferentiable model has been a difficult problem. See Lasdon (1970), Luenberger (1984) and Fletcher (1987) for details about nonlinear programming. Global optimization is a relatively new field which aims at seeking an optimal solution in nonlinear optimization models. The main stream in global optimization exploits analytical properties of the model to generate a deterministic sequence of points converging to a global optimal solution. Two analytical properties, convexity and monotonicity, have been successfully exploited, giving rise to two research trends in global optimization: d.c. optimization and monotonic optimization. The d.c. optimization deals with models described by differences of convex functions, while monotonic optimization deals with models described by means of functions monotonically increasing or decreasing along rays, hence the proposed algorithms still depend on these analytical properties of models. Refer to Horst et al. (1995), Horst and Pardalos (1995), for more about global optimization. Therefore, we can see that only models with some good structures have been solved. For the general nonlinear model, how to find out an optimal solution remains to be a challenging problem. Remark 1. Other methods for solving optimization problems have been proposed that do not necessarily aim at seeking an optimal solution with certainty, such as the stochastic approaches and the metaheuristics. Stochastic approaches generally require little or no additional assumptions on the structure of optimization models, but the outcome they lead to is random, although some provide a probabilistic convergence guarantee, refer to Pardalos et al. (2000) for a review of three kinds of stochastic methods. These approaches are not practically usable since they can not quantify what a result they will produce within a given time period. Metaheuristics are methods that are often based on processes observed in physics or biology. Three main categories of them are: simulated annealing,
20
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
tabu search, and genetic algorithms. We do not dwell on metaheuristics here, but refer readers to Laarhoven and Aarts (1987) for simulated annealing, Glover and Lagunna (1997) for Tabu search, and Goldberg (1989) for genetic algorithms. Many experiments showed that metaheuristics were very successful in solving hard combinatorial optimization problems, but they are heuristic in nature, and they can not quantify the result they lead to. The soft approach to be introduced in this paper differs essentially from the stochastic approaches and the metaheuristics because we can quantify the result obtained, which will be explained in the following sections.
3. A soft approach As stated in the previous section, only a small part of optimization models has been solved. To make matters worse, some ‘‘solved’’ models may not be really solved if taking some practical considerations into account. For instance, when there is a time limit for computation, algorithms which can not lead to an optimal solution within the given time limit will be helpless, because an algorithm converging to an optimal solution in 10 minutes may produce a very bad result when terminated in 8 minutes. Regretfully, most traditional algorithms are not capable to handle time limit requirements. For solving a larger class of optimization models, and also for handling the time limit requirement, this section introduces a soft approach, which does not require any assumption about the objective function, and can handle requirements on computation time. We first explain the ideas in the soft resolution approach, and then suggest a two-stage process for implementing it.
optimization, which was proposed for solving complicated stochastic optimization models, refer to Lau and Ho (1997), and Lee et al. (1999) for details on ordinal optimization. This goal softening method softens the goal of solving an optimization model in two dimensions: 1. Instead of looking for an optimal solution, it seeks a good enough solution. To define a good enough solution, we donÕt use the cardinal performance but the order of a solution in the feasible set. For instance, the top k% solutions may be taken as good enough solutions, where k is set by the requirements in solving a model. 2. Instead of requiring the solution to be an optimal solution for sure, it accepts a solution if it is high likely to be a good enough solution. When the probability with which an alternative is a good enough solution is larger than some value, say 95%, this alternative can be taken as a solution. In other words, the soft approach agrees that a model is ‘‘solved’’ if an alternative which is high likely a good enough solution is obtained. The notions of good enough and high probability are set by the requirements in solving a model. To focus on the main issues, letÕs take the top 1% solutions as good enough, and 95% as high probability. Therefore, when a solution is founded to be in the top 1% with a probability bigger than 95%, the model is solved. The good enough solutions usually form a subset of X f, seeking any point of a set is much easier than seeking exactly one point, it is easy to understand that ‘‘solving’’ an optimization model under the soft approach should be easier. In addition, the soft approach allows failure, i.e., we do not require the solution found to be a good enough solution for sure, it is reasonable to believe that a large class of models are expected to be solvable under the soft approach.
3.1. Ideas in the soft approach 3.2. A two-stage process The soft approach softens the goal in ‘‘solving’’ an optimization model. Among several goal softening methods, we use the one suggested in ordinal
We suggest to implement the soft approach in the following two stages.
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
• Stage 1: Sample the feasible set X f to generate a finite subset S. As for the sample set S, we require the probability with which S contains at least one good enough solution to be high. Let G X f be a set of good enough solutions, say the top 1% solutions, we require the following condition to be satisfied: pfjS \ Gj P 1g P 95%:
ð2Þ
Consequently, the best alternative in S will be a good enough solution with a probability not smaller than 95%. • Stage 2: Select the best alternative x* from the sample set S. Then x* is the solution we are looking for. This two-stage process is illustrated in Fig. 1. Since S is a finite set, there is no conceptual difficulty in Stage 2. However, Stage 2 will be time consuming if S is too big. Fortunately, we will show below that a small sample set is enough to satisfy condition (2), no need to generate a big sample set, so Stage 2 will not take much computation time. Suppose that we generate S by uniformly taking 1000 samples from X f, then the probability that none of the good enough solutions is included in S is p = (1 1%)1000 < 0.1%, consequently, p{jS \ Gj P 1} = 1 p > 99.9%. That is, a sample set with 1000 uniform samples can satisfy (2). In other words, if we generate 1000 uniform samples from X f in Stage 1, and choose the best sample in Stage 2, then the best sample is one of the top 1% solutions with a probability not smaller than 99%, so the model is solved.
21
3.3. Features of the soft approach From the two-stage process, we can easily see that the soft approach has the following two distinguished features. • No assumption about the objective function is needed. The objective function is used only in Stage 2. Since S is a finite and small set, say with 1000 samples, we can evaluate all samples in S to get the best sample, so no assumption about the objective function is needed. • Time for solving a model is controllable. We can estimate the computation time needed for solving an optimization model given the requirements in solving the model, and quantify the quality of the solution obtained within a given time limit. Let t1 be the time for taking one feasible sample, and t2 be the time for evaluating one sample, then the time for solving a model, denoted by T, is equal to T ¼ jSjðt1 þ t2 Þ;
ð3Þ
where jSj is the number of samples taken in Stage 1, which should be determined by the requirements in solving an optimization model and the sampling method used in Stage 1. Suppose it is required that the solution should be one of the top k% alternatives with a probability not smaller than q%, and uniform sampling is used in Stage 1, then jSj can be obtained by solving the following equation: jSj
O O
O O O O
O O O O O
S ** * * * * * ** O * G O
O
G: Set of good enough solutions S: Set of samples O Sampling
O
O Selecting O O O
*
** * S o: Feasible solution *: Good enough solution
Fig. 1. A two-stage process.
pfjS \ Gj P 1g ¼ 1 ð1 k%Þ
¼ q%:
ð4Þ
T is known by substituting jSj into (3). So we can estimate the time needed for solving a model when requirements of solving a model are given. For instance, suppose the solution obtained should be one top 1% alternative with a probability bigger than 99.99%, then we know from (4) that the minimal number of uniform samples we need to take in Stage 1 is logð199:99%Þ ¼ 917. If generatlogð11%Þ ing one uniform sample takes 103 seconds and
22
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
calculating the value of objective function at each sample point takes 104 seconds with a computer, then total computation time needed to solve the problem with this computer is approximately 917 * (103 + 104) 1 second. On the other hand, we can quantify the quality of the solution obtained given a time limit for solving a model. Let T0 be the time limit given, then the number of samples we can take in Stage 1 is jSj ¼
T0 : t1 þ t2
ð5Þ
When Stage 1 takes samples with some uniform sampling method, substituting (5) into the following equation: jSj
pfjS \ Gj P 1g ¼ 1 ð1 k%Þ ;
ð6Þ
we can know the probability with which the solution is one of the top k% solutions. For instance, suppose we are given 2 seconds for solving a model, and it takes a computer 103 seconds for generating one uniform sample and 104 seconds for calculating the value of objective function at each sample point, then the maximum number of samples we can take in Stage 2 1 is 103 þ10 4 ¼ 1818. If the good enough solutions are defined as the top 0.1% points, then the soft approach can generate a good enough solution in 2 seconds with a probability of 1(10.1%)1818 = 83.78%, thus the quality of the solution obtained within 2 seconds is quantified. Remark 2. Time needed for taking one feasible sample depends on three factors: the sampling method in Stage 1, the feasible set of a model and computer speed. Given these three factors, we can know t1 by doing a simple sampling test. Time needed for evaluating one sample depends on the speed of a computer and the objective function of a model, we anticipate that t2 will be much smaller than t1, so time needed in solving an optimization model mainly depends on t1. Computing time controllability may be most valuable feature of the soft approach, cause other optimization methods, including traditional optimization methods, stochastic proaches and metaheuristics, can not handle
the bethe apthe
requirement on computing time, and there are many situations where models need to be solved within some time limit. The soft approach not only gives a solution within a give time limit, but also tells the quality of the solution obtained.
4. An algorithm for solving continuous models This section presents an algorithm for solving continuous models following the two-stage process. To complete Stage 1, we need to decide the number of necessary samples and a method for taking samples from X f. Suppose that the requirement for solving a model is: seeking one of the top k% solutions with a probability not small then q%. Then the sample set S needs to satisfy the following condition: pfjS \ Gj P 1g P q%:
ð7Þ
Suppose that no information is available about the distribution of the top k% solutions, and Stage 1 uses some uniform sampling method. Then the smallest number of necessary samples is determined by solving the following equation with jSj as the variable: 1 ð1 k%ÞjSj ¼ q%:
ð8Þ
We next consider methods for uniform sampling. Let ð9Þ X f ¼ fx 2 X jgj ðxÞ 6 0; j ¼ 1; . . . ; mg; Qn where X ¼ i¼1 X i is the domain of x. Without loosing generality, we assume that points satisfying all constraints are also in the domain.1 There proposed two approaches for uniformly sampling a set described by (9): the rejection approach or the random generating approach. The rejection approach takes samples from X uniformly and accepts samples if they are also in X f and reject otherwise. This method applies to uniformly sample any bounded sets. However, when the feasible set X f is just a very small part of the domain, a lot of points will be rejected be1
This means that X1 = {x 2 Rnjgi(x) 6 0, i = 1,. . ., m} X holds. When X1XX, we can add constraints like ai 6 xi 6 bi to X1 so that X1 X will hold.
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
fore a feasible point is founded, Stage 1 will be time consuming due to the low sampling efficiency. The random generating approach suggests generating uniform samples within the feasible set, its efficiency is high if a proper generating method can be worked out. Smith (1984) proposed a general mixing method for generating uniform samples in a bounded set, which is stated below. • 1. 2. 3.
A general mixing method Find a start point x0 2 X f, let i = 0. Generate a random direction d in Rn. Form a line set L = X f \ {xjx = xi + kd}, where k is a real scalar. 4. Generate a random point xi+1 uniformly distributed over L. 5. Set i to i + 1, go back to Step 2 till necessary samples are obtained.
This general mixing method is illustrated in Fig. 2. This mixing method is good in generality, but it did not tell how to determine k in Step 3. Berbee et al. (1987) proposed a method called ‘‘hit-and-run’’ algorithm for uniformly sampling a simplex. The main difference with the above algorithm is in Step 3, the hit-and-run algorithm forms the line L by connecting two hit points in the boundary of a simplex, which are decided by hitting the boundary of a simplex along directions d and d. We employ this idea to develop an algorithm for uniformly sampling a convex set. The main idea of our algorithm is as follows: Since X f is convex, hitting the boundary of X f from an interior point along one direction will produce only one feasible hit point, and the line seg-
ment connecting any two feasible hit points is still in X f. For a uniformly generated direction d, we can get two hit points in boundary, one in direction d and the other in direction d, then the point uniformly distributed on the line segment connecting the two hit points is a uniform sample. The algorithm is stated in the following. • A uniform sampling algorithm for convex sets 1. Find an interior point x0 2 X f uniformly as a start point, let i = 0. x0 can be founded by using the rejection method while uniformly sampling X. 2. Select a direction d uniformly over a unit sphere D = {d 2 Rnjjjdjj = 1}.2 3. Find the two hit points by hitting the boundary of X f from standing point xi along directions d and d. Denote Pi, Qi: hit points in the boundary of X f along directions d and d, respectively, k1: distance between xi and Pi, k2: distance between xi and Qi. Then we have P i ¼ xi þ k1 d;
d xi
x
i+1
Qi ¼ xi k2 d: 1
ð10Þ
2
We can calculate k and k as follows: Find hit points in each constraint surface gj (x) = 0, wthe distance between xi and the nearest hit point in the constraint surface along direction d is k1, while the distance along direction d is k2. For each j 2 {1,. . .,m}, denote k1j : k2j :
L
23
distance between xi and the nearest hit point in constraint surface gj (x) = 0 along direction d, distance between xi and the nearest hit point in constraint surface gj (x) = 0 along direction d.
Then we know that k1j is the smallest positive solution of the following equation:
Xf Fig. 2. A general mixing process.
2 Taking uniform samples from nice regions like a hyper cube or a hypersphere is easy, refer to Rubinstein (1982).
24
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
gj ðxi þ kdÞ ¼ 0;
ð11Þ
while k2j is the smallest positive solution of gj ðxi kdÞ ¼ 0;
ð12Þ
When (11) or/and (12) do not have a position solution, set k1j or/and k2j to 1, thus we have k1 ¼ minfk11 ; ; k1m g; k2 ¼ minfk21 ; ; k2m g:
For the case where all constraints are linear, finding the hit points is quite easy. Refer to Berbee et al. (1987) for details. When necessary samples are obtained using the above sampling algorithm, Stage 2 selects the best alternative by evaluating all samples in S, leading to one alternative which meets the requirements in solving the model.
ð13Þ
i+1
uniformly over the line segment 4. Choose x connecting the two hit points: xiþ1 ¼ xi þ lðP i Qi Þ;
ð14Þ
where l is a random variable distributed uniformly over [0,1]. 5. Set i to i + 1, go back to Step 2 till necessary samples are obtained. This algorithm is illustrated in Fig. 3. Remark 3. To ensure getting two hit points in each iteration, it is necessary to hit the boundary from an interior point, two strategies for doing this are: (1) Always hit the boundary from the start point x0, (2) If a hit point in the boundary is chosen in Step 4 (l = 0 or 1), do not change standing point in next hitting. When there is not an interior point in X f, refer to Xu and Ho (2002) for a uniform sampling method.
5. Numerical examples We will use the soft approach to solve two optimization problems for illustrating the validity of the soft approach. Example 1. Find the minimum of the following function in interval [2, 2]: ( xj sin ð1xÞj ; x 6¼ 0; x2 þ1 f ðxÞ ¼ ð15Þ 0; x ¼ 0: The figure of function (15) is shown in Fig. 4. This is a sort of hard optimization problem since the objective function is not differentiable at many points, and there are many local minimums. For testing the soft approach, we solved this problem with some other method, the minimum is reached at x0 = 0.7374 and f (x0) = 0.4667. To see how well the soft approach solves this model, suppose that the set of good enough
Remark 4. Seeking hit points is the main computation load in this sampling algorithm, which requires to solve Eqs. (11) and (12) for each constraint function with k as the variable. Since these equations contain only one variable k, there should be no essential difficulty in solving these equations.
0.5 0.4 0.3 0.2
f(x)
0.1 0 -0.1
P
1
d
Q
i
i
x i+1
2
d
-0.2 -0.3 -0.4
xi X
f
Fig. 3. A uniform sampling algorithm for convex sets.
-0.5 -2
-1.5
-1
-0.5
0
0.5
1
1.5
x
Fig. 4. The objective function in Example 1.
2
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
G fx 2 Rjjx x0 j 6 0:04g:
ð16Þ
f
Since X is an interval [2, 2], uniformly sampling X f is very easy, we do not need to use the sampling algorithm proposed in Section 4. After getting 1000 uniform samples in Stage 1, we choose the best sample as the solution. We solved this optimization problem 50 times with the soft approach, and got one good enough point each time. Thus we have illustrated the effectiveness of the soft approach in this example. Table 1 below summarizes the results produced, wherein x* is the solution obtained, and f (x*) is the value of the objective function at point x*. For showing the results clearly, we plot in Fig. 5 the objective values of the best samples in the 50 experiments.
-0.4666 -0.4666 -0.4666 -0.4667
f(x*)
solutions G is defined as the top 1% points in X f = [2, 2], and the solution produced should be a good enough solution with a probability bigger than 99.99%. When taking 1000 uniform samples from X f, the possibility that at least one point in G is included in the samples is p(jG \ X fj P 1) = 10.991000 = 99.996%, so taking 1000 uniform samples in Stage 1 is enough. To check if the solution obtained from using the soft approach is in G, we need to determine G and do some computation experiments. Let the length of an interval be its metric, then the metric of G, denoted by jGj, is 1% * jX fj = 0.04. Fig. 4 shows that the objective function is nearly symmetric in the neighborhood of x0, we can take the points in the neighbor of x0 as an approximation of G, that is,
25
-0.4667 -0.4667 -0.4667 -0.4667 -0.4667 0
5
10
15
20
25
30
35
40
45
Fig. 5. Objective values of the best samples in 50 experiments (Example 1).
We can also examine the degree of approximation in term of value since we know the optimal solution. Let the degree of approximation be measured by f ðx Þ f ðx0 Þ ; d ¼ f ðx0 Þ where x* is a solution obtained from using the soft approach, x0 is the optimal solution. We calculated d in each of the 50 experiments, the worst approximation turned out to be 0.0122%, so the soft approach leaded to a good enough solution even if in term of value. Example 2. Find the minimum of the following function: f ðx1 ; x2 Þ ¼
j sinð5px1 Þ þ sinð5px2 Þj e10ðx1 0:25Þ
2
þ10ðx2 0:5Þ2
ð17Þ
;
Table 1 Computation results from solving Example 1 x* f (x*) x* f (x*) x* f (x*) x* f (x*) x* f (x*)
0.7384 0.4667 0.7305 0.4666 0.7380 0.4667 0.7376 0.4667 0.7369 0.4667
0.7380 0.4667 0.7353 0.4667 0.7374 0.4667 0.7376 0.4667 0.7412 0.4667
0.7370 0.4667 0.7393 0.4667 0.7362 0.4667 0.7375 0.4667 0.7363 0.4667
0.7368 0.4667 0.7364 0.4667 0.7374 0.4667 0.7396 0.4667 0.7409 0.4667
0.7386 0.4667 0.7386 0.4667 0.7380 0.4667 0.7392 0.4667 0.7384 0.4667
50
Sequence of experiments
0.7378 0.4667 0.7384 0.4667 0.7385 0.4667 0.7407 0.4667 0.7361 0.4667
0.7334 0.4667 0.7363 0.4667 0.7384 0.4667 0.7339 0.4667 0.7368 0.4667
0.7324 0.4667 0.7341 0.4667 0.7368 0.4667 0.7302 0.4666 0.7404 0.4667
0.7444 0.4666 0.7380 0.4667 0.7382 0.4667 0.7379 0.4667 0.7387 0.4667
0.7378 0.4667 0.7394 0.4667 0.7402 0.4667 0.7380 0.4667 0.7354 0.4667
26
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
in the feasible set defined by 2
X f ¼ fðx1 ; x2 Þ 2 ½0; 1 j0 6 1 x1 x2 6 1; x1 þ 2x2 P 1g:
ð18Þ
The figure of the objective function (17) is shown in Fig. 6, and X f is shown in Fig. 7. This is a sort of hard problem since the object function is not convex. For testing the soft approach, we solved this problem with some other method, the minimum is reached at (0.293, 0.5) and f (x0) = 1.9574, that is, x0 = (0.293, 0.5), f (x0) = 1.9574. Similar as in Example 1, let the set of good enough solutions be the top 1% points in X f, we
0
f(x1,x2)
-0.5 -1 -1.5
are satisfied with obtaining one good enough solution with a probability bigger than 99.99%, hence we need only taking 1000 uniform samples from X f in Stage 1. Fig. 7 shows that X f accounts for only 1/4 of the domain [0,1]2, we use the sampling algorithm proposed in Section 4 for sampling X f efficiently. Fig. 8 shows the 1000 samples taken uniformly from X f with the proposed algorithm. After getting 1000 uniform samples, choose the best sample as the solution. In order to check if the soft approach produces a good enough solution with a high probability, we need to determine G and do computation experiments. Let the area of a set in R2 be its metric, then jX fj = 13/4 = 0.25, and jGj = jX fj * 1% = 0.0025. Fig. 6 shows that the objective function is nearly symmetric in the neighborhood of x0, we take the circle around x0 as an approximation of G, that is, G {x 2 R2jjjxx0jj 6 r}. Since jGj = 0.0025, r is determined from solving equation pr2 = 0.0025, that is, r = 0.00282, hence we have G fx 2 R2 jjjx x0 jj 6 0:0282g:
-2 1 0.8
1
0.6 0.4
x2
0.2 0
0.2
0
0.4
0.6
0.8
x1
ð19Þ
We solved this optimization problem 50 times with the soft approach, and got one good enough point each time. Thus we have illustrated the effectiveness of the soft approach in this example.
1
1
0.9
0.9
0.8
0.8
0.7
0.7
0.6
0.6
X2
X2
Fig. 6. The objective function in Example 2.
0.5
0.5
0.4
0.4
0.3
0.3
0.2
0.2
0.1
0.1
0
0
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
X1
Fig. 7. Feasible set in Example 2.
0.8
0.9
1
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
X1
Fig. 8. Uniform samples from X f.
0.8
0.9
1
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
27
Table 2 Computation results from solving Example 2 x1 x2 f (x*) x1 x2 f (x*) x1 x2 f (x*) x1 x2 f (x*) x1 x2 f (x*)
0.2886 0.4997 1.9546 0.2949 0.5008 1.9568 0.2882 0.4960 1.9519 0.3009 0.5031 1.9475 0.2898 0.4920 1.9471
0.2958 0.5084 1.9465 0.2861 0.4982 1.9501 0.2928 0.5062 1.9519 0.2978 0.5018 1.9537 0.2987 0.5086 1.9424
0.2867 0.4913 1.9413 0.2949 0.4909 1.9454 0.2979 0.4956 1.9514 0.2878 0.5045 1.9507 0.2917 0.4975 1.9563
0.2993 0.5045 1.9490 0.2897 0.4962 1.9539 0.2881 0.4927 1.9465 0.2866 0.4921 1.9429 0.2946 0.4965 1.9554
0.3013 0.5042 1.9455 0.2873 0.5016 1.9524 0.3029 0.4956 1.9411 0.2915 0.5061 1.9519 0.2874 0.4996 1.9530
Table 2 summarizes the results produced, wherein x* is the solution obtained, and f (x*) is the value of the objective function at point x*. For seeing the results clearly, we plot in Fig. 9 the objective values of the best samples in the 50 experiments. We can also examine the degree of approximation in term of value of since we know the optimal solution. We calculated the degree of approximation in each of the 50 experiments, the worst degree of approximation turned out to be 1.9752%, so the soft approach leaded to a good solution even if in term of value. -1.915 -1.92 -1.925 -1.93
f(x*)
-1.935 -1.94 -1.945 -1.95 -1.955 -1.96
0
5
10
15
20
25
30
35
40
45
50
Sequence of experiments
Fig. 9. Objective values of the best samples in 50 experiments (Example 2).
0.2933 0.5168 1.9182 0.2895 0.5013 1.9554 0.2885 0.4981 1.9541 0.2870 0.5003 1.9523 0.2898 0.4991 1.9559
0.3096 0.5005 1.9190 0.3015 0.5039 1.9451 0.2978 0.5059 1.9494 0.2946 0.5012 1.9569 0.2815 0.4986 1.9384
0.2854 0.5029 1.9481 0.2979 0.5001 1.9541 0.2907 0.5019 1.9561 0.3057 0.4954 1.9322 0.3062 0.5049 1.9299
0.2831 0.5008 1.9435 0.2954 0.5052 1.9529 0.2875 0.4920 1.9442 0.2815 0.4907 1.9263 0.2960 0.5008 1.9561
0.2973 0.5060 1.9498 0.2864 0.5022 1.9506 0.2982 0.5028 1.9526 0.2997 0.4874 1.9292 0.2933 0.4966 1.9558
Remark 5. We specified the good enough set in the two examples, this is only for checking if the result obtained is G. The soft approach does not require that we should specify G, and it is also not necessary to know where the good enough solutions are located in the feasible set in using the soft approach. We got the optimal solution with other optimization methods in the two examples, this is for specifying the good enough set. It is not necessary to know the optimal solution for quantifying the solution produced by the soft approach, because the soft approach does not use value in quantifying the solution produced. Computation time of the soft approach is determined by the time needed for generating necessary samples and calculating the values of the objective function at the sample points, which depends on the speed of computer used. We used a IBM ThinkPad computer (CPU: Pentium M Processor 1.6 GHz, Memory: 768MB) in solving the above model. It takes 1.5 * 104 seconds in generating one uniform sample with the proposed method, and 7 * 105 seconds in evaluating one sample, so solving this model takes 1000 * (1.5 * 104 + 7 * 105) = 0.22 seconds. To see what a solution the soft approach can produce within a time limit, we do some computing experiments by limiting the computing time to 2 seconds.
28
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
With the same computer, the maximum number samples we can take in Stage 1 is 2 ¼ 9091. Suppose the good enough 1:5104 þ7105 solutions are defined as the top 0.1% points, then the soft approach is expected to produce within 2 seconds one good enough solution with a probability of 1(10.1%)9091 = 99.99%. To verify this theoretical anticipation, we specify G first and then do 50 computing experiments. Since jGj = 0.1%jX fj = 0.00025 and G {x 2 R2jjjxx0jj 6 r}, we can know r from solving equation pr2 = 0.00025, that is, r = 0.0089, hence of
G fx 2 R2 jjjx x0 jj 6 0:0089g:
ð20Þ
We solve this optimization problem by taking 9091 uniform samples in Stage 1 and choosing the best sample as the solution in Stage 2. We repeated this solving process 50 times, and got one point in G defined by (20) each time, thus the soft approach did produce a good enough solution within a given time limit. The objective values of the best samples obtained from the 50 experiments are shown in Fig. 10. We omit the detailed experiment data for saving space. We calculated the degree of approximation in each of the 50 experiments, and the worst degree of approximation turned out to be 0.31%. Comparing with the worst degree of approximation
6. Conclusions This paper introduced a soft approach for solving complicated optimization models, and suggested a two-stage process for implementing the soft approach. An algorithm was given for solving optimization problems with a convex feasible set. The soft approach has two distinguished features that most traditional optimization methods do not have, i.e., no assumption about the objective function is needed, and the computation time is controllable. The soft approach differs essentially from stochastic and heuristic methods in the sense that the result obtained can be quantified. The soft approach does not intend to replace the optimum seeking methods, it does provide an alternative when optimum seeking fails. We expect that the soft approach will be helpful when no method is available for seeking an optimal solution, or when existing algorithms can not generate a good solution within a given time limit. Further research on the soft approach can be conducted in the following directions: (1) Efficient sampling methods. Sampling plays the central role in the soft approach, designing efficient sampling methods for different sets is the key for getting a better result. (2) Biased sampling methods. When there is some information about the distribution of good enough solutions, concentrating on the areas with a high density of good enough solutions in sampling will produce better results than uniform sampling, how to design such sampling methods is an interesting issue.
-1.951
-1.952
-1.953
-1.954
f(x*)
when taking 1000 samples, which was 1.9752% in our previous 50 experiments, we can see that taking more samples in Stage 1 improved the quality of solution even in term of value.
-1.955
-1.956
-1.957
-1.958
Acknowledgments 0
5
10
15
20
25
30
35
40
45
50
Sequence of experiments
Fig. 10. Objective values of the best samples in 50 experiments when taking 9091 samples.
Dr. C. Xu is grateful to Professor Y.C. Ho at Harvard University for his kind support and insightful guidance while doing research with him
C. Xu, P. Ng / European Journal of Operational Research 173 (2006) 18–29
at Harvard, and also for the wonderful ideas in his ordinal optimization, from which this study benefits a lot.
References Berbee, H.C.P., Boender, C.G.E., Kan, A.H.G.R., Scheffer, R.L., Smith, R.L., Telgen, J., 1987. Hit-and-run algorithms for the identification of nonredundant linear inequalities. Mathematical Programming 37 (2), 184–207. Danzig, G.B., 1963. Linear Programming and Extensions. Princeton University Press, Princeton, NJ. Fletcher, R., 1987. Practical Methods of Optimization. John Wiley & Sons, New York. Glover, F., Lagunna, M., 1997. Tabu Search. Kluwer Academic Publishers, Boston, MA. Goldberg, D., 1989. Genetic Algorithms in Search, Optimization, and Machine Learning. Addision-Wesley, San Mateo, CA. Horst, R., Pardalos, P.M., Thoai, N.V., 1995. Introduction to Global Optimization. Kluwer Academic Publishers, Dordrecht. Horst, R., Pardalos, P.M. (Eds.), 1995. Handbook of Global Optimization. Kluwer Academic Publishers, Dordrecht.
29
Laarhoven, P.J.M. van, Aarts, E.H., 1987. Simulated Annealing: Theory and Applications. Kluwer Academic Publishers, Dordrecht. Lasdon, L.S., 1970. Optimization Theory for Large Systems. Macmillan, London. Lau, T.W.E., Ho, Y.C., 1997. Universal alignment probabilities and subset selection for ordinal optimization. JOTA 93 (3), 455–489. Lee, L.H., Lau, T.W.E., Ho, Y.C., 1999. Explanation of goal softening in ordinal optimization. IEEE Transactions on AC 44 (1), 94–99. Luenberger, D.G., 1984. Introduction to Linear and Nonlinear Programming. Addison-Wesley, Reading, MA. Pardalos, P.M., Romeijn, H.E., Tuy, H., 2000. Recent developments and trends in global optimization. Journal of Computational and Applied Mathematics 124, 209–228. Rubinstein, R.Y., 1982. Generating random vectors uniformly distributed inside and on the surface of different regions. European Journal of Operational Research 10 (2), 205–209. Smith, R.L., 1984. Efficient Monte-Carlo procedures for generating points uniformly distributed over bounded regions. Operations Research 32 (6), 1296–1308. Xu, C., Ho, Y.C., 2002. Constrained ordinal optimization, Unpublished research report, Division of Engineering and Applied Sciences, Harvard University.