INFLUENCE OF LIPSCHITZ BOUNDS ON THE SPEED OF GLOBAL OPTIMIZATION

Global optimization methods based on Lipschitz bounds have been analyzed and applied widely to solve various optimization problems. In this paper a bound for Lipschitz function is proposed, which is computed using function values at the vertices of a simplex and the radius of the circumscribed sphere. The efficiency of a branch and bound algorithm with proposed bound and combinations of bounds is evaluated experimentally while solving a number of multidimensional test problems for global optimization. The influence of different bounds on the performance of a branch and bound algorithm has been investigated.


Introduction
Global optimization is used in such fields as economics, operations research, biology, engineering design and management as well as other numerous engineering and applied sciences. For example, prediction methods in computational finance forecast various kinds of quantities related to stock markets, like stock prices, stock volatility and ranking measures. Prediction problem of building a multi-stock artificial trader (Björkman, Holmström 1999) can be used instead of the classical time series approach. The behavior of the trader is controlled by a parameter vector which is tuned for the best performance. Global optimization algorithm DIRECT (Jones et al. 1993) was applied to find optimal parameters. Extension of optimization to decision-making systems is a challenging topic of research (Sakalauskas, Zavadskas 2009).
One of the most investigated fields of global optimization is Lipschitz optimization (Strongin 1978;Strongin, Sergeyev 2000;Sergeyev, Kvasov 2008). In this paper we consider the problem of global optimization of a Lipschitz continuous objective function : → x differs from * f by no more than a specified accuracy ε can be found. There are at least four approaches to obtain an estimate of Lipschitz constant (Sergeyev, Kvasov 2009): 1) it can be given a priori; 2) its adaptive global estimate over the whole domain can be used; 3) local Lipschitz constants can be estimated adaptively using the local tuning technique (Sergeyev 1995(Sergeyev , 1998; 4) a set of possible estimates can be used (Jones et al. 1993).
The most studied case of problem (2) is the univariate one ( =1) n , for which numerous algorithms have been proposed, compared, and theoretically investigated. In the present paper, we are interested in the multivariate case ≥ ( 2) n . In Lipschitz optimization the upper bound for ( ) f x over a sub-region ⊆ I D is evaluated by exploiting Lipschitz condition (1). It follows from (1) that for all ∈ , If ∈ y I is fixed, then the convex function overestimates ( ) f x over I . The superscript q indicates the norm used for the upper bound calculation. Let T be a finite set of distinct points in I . Then, the Piyavskii-type upper bound over I , given the function values ∈ ( ), f y y T , and the Lipschitz constant L , is provided by In the univariate case, the function ∈ ( ) min q y y T F x is piecewise linear, and ϕ q can be determined in a simple straightforward way (Piyavskii 1972). For (n > 2), however, problem (4) constitutes a difficult optimization problem.
Therefore, branch and bound algorithms use considerably weaker bounds. In general, weaker bounds belong to the following two simple families µ 1 and µ 2 . Let x y x y denote the diameter of I . For example, if = { : , and if I is an n -simplex, then the diameter δ ( ) q I is the length of its longest edge. Afterwards a simple upper bound can be derived from (3): where ⊂ T I is a finite sample of points in I , where the function values of f have been evaluated. If I is a rectangle or a simplex, the set T often coincides with the vertex set ( ) V I . A tighter but computationally more expensive than (5) bound is Methods with the ϕ q type bounds (4) can hardly be used to solve typical test problems with relatively large Lipschitz constants and (n > 2). Branch and bound algorithms can almost always provide reasonable approximate optimal solutions for (n = 3) , although the methods usually involve more function evaluations (and thus are less suitable in case of very expensive functions), but much less auxiliary computational time than the methods with the ϕ q type bounds. Therefore, it is important to investigate possibilities of bounds tighter than (5) and (6), but computationally less expensive than (4). In this paper, we propose a bound based on the function values at the vertices and the circumradius of an n-circumsphere.

Branch and bound with simplicial partitions and Lipschitz bound
Branch and bound is a technique for implementation of covering global optimization methods as well as combinatorial optimization algorithms. An iteration of a classical branch and bound algorithm processes a node in the search tree representing a not yet explored subspace of the solution space. The iteration has three main components: selection of the node to process, branching of the search tree and bound calculation. Several selection strategies (Paulavičius et al. 2010) may be used: best-first, depth-first, breadth-first, statistical (Žilinskas, A., Žilinskas, J. 2010). Partitions may be hyper-rectangular, simplicial (Horst 1997;Gorodeckij 1999), hyper-conic or hyper-spherical. Bounds may be computed using envelopes of functions (Androulakis et al. 1995;Adjiman, Floudas 1996;Zlobec 2010), interval arithmetic (Moore 1966;Hansen, Walster 2003;Žilinskas 2005) as well as Lipschitz condition. Most of Lipschitz global optimization branch and bound algorithms use hyper-rectangular partitions. Simplicial partitions are preferable when the values of an objective function at the vertices of partitions are used to compute bounds. Another advantage of simplicial partitions is that they may be used to vertex-triangulate feasible regions of non-rectangular shape defined by linear inequality constraints (Žilinskas 2008), what allows the reduction of search space of problems with symmetric objective functions (Žilinskas 2007).
For simplicial branch and bound, the feasible region should be initially covered by simplices. The most preferable initial covering is face to face vertex triangulation -partitioning of the feasible region into finitely many n-dimensional simplices, whose vertices are also the vertices of the feasible region. We use a standard way of triangulation into n! simplices. All simplices share the diagonal of the feasible region and are of equal hyper-volume. It is also possible to over-cover the feasible region by one simplex (Žilinskas, A., Žilinskas, J. 2002). However, the objective function may be undefined outside the feasible region. Moreover, the hyper-volume of over-covering simplex is considerably larger than the feasible region. For example, one version of over-covering is to fit a hyper-rectangle into a simplex matching a vertex: one vertex of the hyper-rectangle and one vertex of the simplex are matched, the edges of the simplex from this vertex include edges of the hyper-rectangle from this vertex, and the opposite vertex of the hyper-rectangle is placed in the opposite face of the simplex. If a hyper-cube is over-covered using this strategy, the hyper-volume of the simplex is exponentially larger than the hyper-volume of the hyper-cube (approximately π / 2 n e n times). There are some possibilities to partition a simplex into smaller ones (Žilinskas, A., Žilinskas, J. 2002). Experiments have shown that the most preferable partitioning is a subdivision of the simplex into two by a hyper-plane passing through the middle point of the longest edge and the vertices which do not belong to the longest edge. The lower and upper bounds for the maximum of the function over a simplex are estimated using the function values at the vertices. At the beginning of the algorithm the rectangular feasible region is covered by simplices and some initial values are initialized. Then the loop is executed until the candidate list is not empty. The simplex with the best upper bounding function value is chosen. If the difference between the upper bound and the global lower bound is larger than the predefined precision ε the simplex is subdivided into two new simplices which are inserted into candidate list. Convergence of such a branch and bound algorithm for maximizing a Lipschitz functions over n-simplices by using bisection at the midpoint of one of the longest edges and the upper bounds µ 1 (5) q , µ 2 (6) q follows from exhaustiveness of the subdivision and continuity of the function f(x) .

Lipschitz bound over simplices based on the function values at the vertices and the radius of the circumscribed sphere
In this section a new Lipschitz bound over simplices is proposed. The bound is often stronger than usually used trivial bounds and still computationally cheap, especially for low-dimensional problems ≤ ( 3) n for which calculation of the radius of the circumscribed sphere is cheap. To find the radius, n + 2 determinants of (n + 1) ×(n + 1)-dimensionality matrices must be calculated.
be a Lipschitz continuous objective function, ⊂ I D is an n-simplex, 2 ( ) R I denotes the radius of circumscribed n-sphere and ( ) V I denotes the set of vertices. Then Proof. An n-simplex I is covered by n-balls Ο + , =1, , 1  i i n such that the radius 2 ( ) R I is the same for all n-balls and centers coincide with the vertices i v of the simplex I (two dimensional example is shown in Fig. 1 Fig. 1. Covering of a triangle (two-dimensional simplex) I by disks (2-balls) Ο + , =1, , 1  i i n of the same radius 2 ( ) R I
Completing the square for (9) gives: which is the n-circumsphere

Combination of bounds with various norms
Experiments have shown (Paulavičius, Žilinskas 2006(Paulavičius, Žilinskas , 2007 that no single norm and corresponding Lipschitz constant is the best for all problems and therefore a combination of the infinite, Euclidean and first norms is used for calculation of the upper bound The Piyavskii type (4) Lipschitz bound with the first norm ϕ 1 was proposed in (Paulavičius, Žilinskas 2008): However, the first norm does not always give the best bounds (Paulavičius, Žilinskas 2008). The further investigation (Paulavičius, Žilinskas 2009) has shown that an aggregate bound (AB) composed of ∞ µ 2, 2 and ϕ 1 bounds

Experimental investigation
In this section the results of computational experiments are presented and discussed. Various test problems ≥ ( 2) n for global optimization from (Hansen, Jaumard 1995;Jansson, Knüppel 1992;Madsen, Žilinskas 2000) have been used in our experiments. The list of all used test problems and their global minima from the literature are shown in Table 1. For the (n = 2, 3) test problems we use the same precision ε as used by Hansen and Jaumard (Hansen, Jaumard 1995). The proposed algorithm assumes that Lipschitz constants are given a priori. These constants have been evaluated in the same way, as in (Hansen, Jaumard 1995). A very fine grid search algorithm of 1000 n points for (n = 2, 3) and 100 n points for (n = 4, 5, 6) dimensional test functions have been used and the obtained estimates should be close to the best Lipschitz constants. The values are given in (Paulavičius, Žilinskas 2006(Paulavičius, Žilinskas , 2007. The speed of global optimization has been estimated using the number of function evaluations criterion.
The results of the simplicial branch and bound algorithm with the proposed bound ψ 2 (7) based on the Euclidean norm and with the bound µ 2 2 (5) based on the same norm are shown in Table 2. For all test problems the number of function evaluations is smaller when the proposed bound ψ 2 is used and on the average it is by 26% smaller than when the bound µ 2 2 is used. The rate ( ψ < µ 2 2 2 ) shows how often the proposed bound ψ 2 is better (more tight) than the bound µ 2 2 . For the used test functions it is better on average in 70% of simplices. However, the bound ψ 2 does not always give better bound comparing with µ 2 2 (the rate is not equal to 1). If function values at the vertices V(I) of the simplex I differ significantly: I , then µ 2 2 type bound is more tight comparing with the proposed ψ 2 bound. If the difference between the maximal and minimal functions values is smaller than ( ) δ − 2 2 2 ( ) ( ) L R I I , then the proposed ψ 2 type bound gives better results. Therefore, for computationally expensive functions it might be worthwhile to append the bound ψ 2 to the aggregate bound AB (13). Performance of simplicial branch and bound algorithm with the aggregate bound AB (13) is similar to that of the best branch-and-bound algorithm for Lipschitz optimization and it is often better (Paulavičius, Žilinskas 2009). The numbers of function evaluations on the average is by 23% smaller when the improved aggregate bound IAB (14) is used than when the aggregate bound AB (13) is used.
The value of objective function at the same point is required when computing bounds over neighbor sub-regions. If this value is evaluated for parent sub-region, it is not necessary to evaluate it again. However, in multidimensional case the same points may be involved in subdivision of different sub-regions. It is possible to maintain a list of all points where the objective function is evaluated and avoid evaluation of objective function at such points several times by verifying if the point is in the list before evaluation. In two-dimensional case after initial covering by simplices there are two simplices. When one of them is subdivided into two by a hyper-plane passing through the middle point of the longest edge the function value at this point should be evaluated. Sometime later in the optimization process subdivision of the other simplex may be performed, again requiring the function value at the same point. By verification of this point it is possible to avoid evaluation of the function again. Such vertex verification reduces the number of function evaluations essentially, but also increases 6,6 0 4L 2 computational time and therefore it is more suitable in the case of computationally expensive functions. The numbers of function evaluations of the algorithm with improved aggregate bound and verification of vertices are shown in Table 2 in the column IAB . These results using improved aggregate bound IAB significantly improve performance of the proposed Lipschitz optimization algorithm (Paulavičius, Žilinskas 2009), whose results were often better than that of the best branch and bound algorithm for Lipschitz optimization. It is interesting to compare the results with other type of bounds. αBB method (Androulakis et al. 1995) assumes that objective function is twice continuously differentiable and uses this feature to construct convex underestimators. Since the assumption on Lipschitz continuity is weaker, in general Lipschitz bounds should be worse than that of αBB, however, the difference is not known. We investigate IAB and a piecewise application of the αBB underestimator (Gounaris, Floudas 2008) to see the difference. Table 3 presents obtained upper bounds for the maximum value of test problem No. 2 after the same number of evaluations of objective function as results given in (Gounaris, Floudas 2008). Their method subdivides the feasible region into N1 × N2 rectangles which requires (N1 + 1) × (N2 + 1) function evaluations. Therefore, we investigate the upper bounds estimated after 4 (N1 = N2 = 1), 9 (N1 = N2 = 2), 25 (N1 = N2 = 3) and so on function evaluations. As expected, αBB overestimator provides better (more tight) upper bounds. In the beginning the difference is larger, but it becomes smaller when the numbers of partitions increase.

Conclusions
In this paper a new Lipschitz bound over simplices is proposed. It is based on the function values at the vertices and the radius of the circumscribed sphere of the simplex. The speed of simplicial branch and bound algorithm with the proposed Lipschitz bounds has been investigated and compared. Test problems of various dimensionalities (n = 2, 3, 4, 5, 6) from the literature have been used for experimental investigation of the algorithm. The experiments showed that the proposed bound is often better than a usually used bound. For the used test functions the rate ( ψ < µ 2 2 2 ) showed that proposed bound ψ 2 is better (more tight) than the bound µ 2 2 on average in 70% of simplices. Therefore, the use of the proposed bound in simplicial branch and bound reduces the number of function evaluations by 26% comparing to a usual bound. If the proposed bound is used in an aggregate bound with various norms performance improvement is 23% comparing to the aggregate bound without the proposed bound.
The proposed method appears to be inferior to the αBB method. This is to be expected because the latter method assumes that the objective function is twice continuously differentiable.