- Email: [email protected]
Applications of interval arithmetic in non-smooth global optimization
Applied Mathematics and Computation 144 (2003) 413–431 www.elsevier.com/locate/amc
Applications of interval arithmetic in non-smooth global optimization Peiping Shen
a,b,*
, Kecun Zhang a, Yanjun Wang
a
a
b
Faculty of Science, Xi’an Jiaotong University, Xi’an, Shaanxi 710049, PR China Department of Mathematics, Henan Normal University, Xinxiang 453002, PR China
Abstract In this paper, an expansion scheme for general functions is presented, which can be used to achieve better enclosure of function ranges. In the context of interval branchand-bound methods for non-smooth global optimization, a pruning technique is given via the expansion. An interval pruning test is established. This pruning test offers the possibility to cut away a large part of the currently investigated box by the optimization algorithm, which can be utilized as an accelerating device similar to the monotonicity test frequently used in interval methods for smooth problems. Numerical computation shows that the proposed method is reliable and can find all solutions of global optimization problem. Ó 2002 Elsevier Inc. All rights reserved. Keywords: Global optimization; Non-smooth function; Branch-and-bound; Interval arithmetic
1. Introduction Theory and methods of global optimization can be found many applications in practice, such as the primary design of optical thin-film, the design of the electron-filter and steam power plant etc. see [11], but due to the global requirement it is not easy in comparison with theory and methods of local optimization. Horst, Tuy and Pardalos have made a progress in this field [1,2]. In these methods for global optimization, some need to compute the gradient of
*
Corresponding author. E-mail address: [email protected] (P. Shen).
0096-3003/02/$ - see front matter Ó 2002 Elsevier Inc. All rights reserved. doi:10.1016/S0096-3003(02)00417-4
414
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
function [3,4], some require only value of function [8,10,12]. They are much more welcomed by practical workers because the given algorithms only need the value of function. This paper presents a method which only compute value of function, we concentrate on the interval branch-and-bound methods (IBBM) for global optimization. IBBM is used to find the solutions of the following global optimization problems min f ðxÞ; x2X 0
f : Rn ! R1
ð1Þ
where X 0 ¼ fx 2 Rn ja0i 6 xi 6 b0i ; ði ¼ 1; . . . ; nÞg is a box and f is a continuous function over X 0 . IBBM often incorporates a monotonicity test (see [3,7]) to delete subboxes from X 0 which do not contain any global minimum point. The test uses an interval evaluation of the derivative or slope of the objective function over the current subboxes. Moreover, the performance of an IBBM is influenced by the magnitude of the interval vector S, where S satisfies expression f ðxÞ 2 f ð~xÞ þ S ðX ~xÞ for any x; ~x 2 X X 0 ;
ð2Þ
which is usually calculated as a gradient or by means of slope. In this paper, we propose rules to evaluate the expansions of non-linear functions, present a new method for computing certified enclosures of the global minimum value and all global minimum points of non-smooth objective functions subject to a box constraint by means of the expansion. The method employs a special pruning test to cut away a large part of the living subboxes, independent of whether the slope interval contains zero or not. This work has the following three advantages in comparison with existing IBBM: (i) The proposed expansion scheme can yield sharper interval vector S and the value of f over a domain than existing methods (see Table 1). (ii) The new pruning technique uses interval expansion scheme without concerning the smoothness of the objective functions. (iii) We do not need the assumption that the objective function has only finitely many local minima, which is made in [9]. This is due to the fact that the pruning test does not delete any global minimizer. Contrary, in [9] it is only proven that the function values in the deleted parts of the current interval are not smaller than the global minimum value. Moreover, the objective function in [9] is only 1-dimensional and at least three times differentiable. So the proposed method is applicable to rather general cases. The organization of the paper is as follows. In Section 2, we give expansion schemes for general functions, and the quality of enclosure for function ranges is improved through various interval intersections and special treatment for locally convex or concave functions. In Section 3, we give an interval pruning test that is used to discard regions in which there exists no global solutions and propose algorithm for finding all solutions of global optimization problem (1). Illustrative examples are presented in Section 4. The paper concludes in Section 5.
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
415
The notations of this paper are as follows: Let IðR1 Þ be the set of real intervals in R1 , for any interval Xi ¼ ½xi ; xi R1 , then IðXi Þ ¼ fYi 2 IðR1 ÞjYi Xi R1 g. The width wðXi Þ and midpoint mðXi Þ of Xi 2 IðR1 Þ are defined by wðXi Þ ¼ xi xi and mðXi Þ ¼ ðxi þ xi Þ=2. For a real interval vector (a box) X ¼ ðXi Þn1 2 IðRn Þ, where Xi ; i ¼ 1; 2; . . . ; n, are real intervals, and the width and midpoint of X are defined by wðX Þ ¼ max1 6 i 6 n wðXi Þ; mðX Þ ¼ ðmðXi ÞÞn1 2 Rn , respectively. The function F : IðRn Þ ! IðR1 Þ is called an interval extension of a function f : Rn ! R1 in X Rn , if x 2 X implies f ðxÞ 2 F ðX Þ. In other words, f ðX Þ ¼ frag ðX Þ F ðX Þ, where f ðX Þ ¼ frag ðX Þ is the range of the function f over X , that is f ðX Þ ¼ ½minx2X f ðxÞ; maxx2X f ðxÞ. It is assumed in the following that the extension functions have the isotonicity property, i.e. X Y implies F ðX Þ F ðY Þ. The global minimum value of f on X 0 is denoted by f , and the set of global minimizer points of f on X 0 by X . Our aim is to propose an algorithm to find X and f in this paper.
2. Expansion of non-linear functions In this section we will establish the rules of expansion arithmetic to yield sharper enclosure for S in (2) and sharper interval extension of f over a given domain. First, we give definitions about expansions of general non-linear functions, then two theorems which show the calculations of S and interval extension of f are presented. e ; X D be Definition 1. Let the set valued function F : D R1 ! IðR1 Þ and X 1 1 1 given. The triple ðFc ; Fr ; Fs Þ 2 IðR Þ IðR Þ IðR Þ is said to expand F in X e , if with respect to X e : 8~x 2 X
F ð~xÞ Fc ;
8x 2 X :
F ðxÞ Fr ;
e ; 8x 2 X ; 8~ 8~x 2 X y 2 F ð~xÞ; 8y 2 F ðxÞ; 9 Fes 2 Fs :
y y~ ¼ Fes ðx ~xÞ:
e is defined by Furthermore, the slope of F with respect to X and X ( ) ~ y y e Þ; x 2 X ; ~x 2 X e ; x 6¼ ~x : e ; X Þ :¼ slðF ; X y 2 F ðX Þ; y~ 2 F ð X x ~x e Þ Fc , F ðX Þ Fr , slðF ; X e ; X Þ Fs , and From Definition 1, we have F ð X e Þ; F ðX Þ; slðF ; X e ; X ÞÞ expands F in X with respect to X e: ðF ð X Definition 2. Let f : X 0 Rn ! R1 be the objective function of problem (1), and e ¼ ðX ei Þ 2 IðX 0 Þ, X ¼ ðXi Þ 2 IðX 0 Þ. For fixed X1 ; . . . ; Xi 1 ; X eiþ1 ; . . . ; X en , X n1 n1
416
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
eiþ1 ; . . . ; X en Þ is a the function f i : R1 ! IðR1 Þ with f i ðyÞ :¼ f ðX1 ; . . . ; Xi 1 ; y; X set valued function of the variable y 2 Xi and is called the component function for f . From Definition 2, we have ei Þ ¼ f ðX1 ; . . . ; Xi 1 ; X ei ; X eiþ1 ; . . . ; X en Þ; f iðX eiþ1 ; . . . ; X en Þ f i ðXi Þ ¼ f ðX1 ; . . . ; Xi 1 ; Xi ; X ei Þ; f i ðXi Þ; slðf i ; X ei ; Xi ÞÞ expands f i in Xi with respect to X ei . and ðf i ð X e For any x 2 X ; ~x 2 X , and any 0 6 i 6 n, it can be inductively obtained that i X ej ; Xj Þ ðXj X ej Þ: ð3Þ f ðx1 ; . . . ; xi ; ~xiþ1 ; ; ~xn Þ 2 f i ðXi Þ f ð~xÞ þ slðf j ; X j¼1
ei Þ ¼ f i 1 ðXi 1 Þ; i ¼ 2; . . . ; n; From Definitions 1 and 2, we have: (i) f i ð X e and (ii) for any fixed ~x 2 X , the right hand side of (3) with i ¼ n is an interval e ; X ÞÞ extension of f over X , and ðslðf i ; X n1 can be viewed as one of S in (2). In the following, we give the rule of expansion arithmetic for the evaluations of component functions in order to compute (3), for i ¼ 1; 2; . . . ; n, let e ¼ ðX ei Þ , X ¼ ðXi Þ IðRn Þ. X n1 n1 If f ðxÞ ¼ kða constant for 8x 2 X Þ, then the triple ðk; k; 0Þ expands its ei , and f i ðXi Þ ¼ k for all component function f i in Xi with respect to X i ¼ 1; . . . ; n. If f ðxÞ ¼ xj ðj ¼ 1; . . . ; nÞ for all x 2 X , then ðfci ; fri ; fsi Þ expands its compoei , where nent function f i in Xi with respect to X ( ( ej if i 6 j; ej if i < j; X X i i e i i fr ¼ f ðXi Þ ¼ and fc ¼ f ð X i Þ ¼ Xj if i > j; Xj if i P j 1 if i ¼ j; e i ; Xi Þ ¼ fsi ¼ slðf i ; X 0 if i 6¼ j: The interval extension for general functions can be constructed by recursively using the results in the following theorems. Theorems 1 and 2 show how the interval extension of a general function can be constructed from known interval extensions of simple functions. In the theorem, the interval extensions of the functions f and g are known, we give the process to construct the interval extension of the function h which is a composite function of f and g, or f and u. e D; X D. f i and gi are the comTheorem 1. Let f ; g : D Rn ! R1 and X ei Þ, f i :¼ f i ðXi Þ, ponent functions of f and g, respectively. Let fci :¼ f i ð X r ei ; Xi Þ, gi :¼ gi ð X ei Þ, gi :¼ gi ðXi Þ, gi :¼ slðgi ; X ei ; Xi Þ. If ðf i ; f i ; f i Þ and fsi ¼ slðf i ; X c r s c r s ðgci ; gri ; gsi Þ expand the component functions f i of f and gi of g in Xi with respect to ei , respectively, then ðhi ; hi ; hi Þ expands the component function hi of the function X c r s ei , where h ¼ f g for 2 fþ; ; ; =g in Xi with respect to X
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
h1c :¼ fc1 gc1 ;
hic :¼ hi 1 r
hir :¼ ðfri gri Þ \
h1c þ
i X
with i 6¼ 1
417
for 2 fþ; ; ; =g; !
ej Þ hjs ðXj X
for 2 fþ; ; ; =g;
j¼1
8 i f gi > > < s i si ðfs gr þ fri 1 gsi Þ \ ðfri gsi þ fsi gri 1 Þ; his :¼ > with gr0 :¼ gc1 and fr0 :¼ fc1 > : i ðfs hic gsi Þ=gri
for 2 fþ; g; for ¼ ; for ¼ =;
provided no division by zero occurs. If uðfcn [ frn Þ D with u : D R1 ! R1 , then for the composite function hðxÞ ¼ uðf ðxÞÞ, we have the component functions h1c :¼ uðfc1 Þ;
hic :¼ hri 1 for i 6¼ 1;
hir :¼ uðfri Þ \
h1c þ
i X
his :¼ slðu; fci ; fsi Þ fsi ; !
ej Þ : hjs ðXj X
j¼1
Proof. We only give the proof for f g, f =g, and uðf Þ, the results for other cases can be similarly proved. e1 Þ ¼ f 1 ð X e1 Þ g1 ð X e1 Þ f 1 g1 and For hðxÞ ¼ f ðxÞ gðxÞ, we have h1 ð X c c i e i 1 1 1 1 i i 1 h ð X i Þ ¼ h ðXi 1 Þ for i 6¼ 1, and hence hc ¼ fc gc , hc ¼ hr . Moreover, ( ) f i ðxi Þ gi ðxi Þ f i ð~xi Þ gi ð~xi Þ i e ei ; xi 2 Xi ; ~xi 6¼ xi : slðh ; X i ; Xi Þ ¼ ~xi 2 X xi ~xi ei ; xi 2 Xi , there exist f~i 2 f i , g~i 2 gi with For any fixed ~xi 2 X s s s s f i ðxi Þgi ðxi Þ f i ð~xi Þgi ð~xi Þ ¼ ½f i ð~xi Þ þ f~si ðxi ~xi Þ½gi ð~xi Þ þ g~si ðxi ~xi Þ f i ð~xi Þgi ð~xi Þ ¼ ½f~i gi ðxi Þ þ f i ð~xi Þ g~i ðxi ~xi Þ: s
s
Similarly, we have f ðxi Þg ðxi Þ f ð~xi Þgi ð~xi Þ ¼ ½f i ðxi Þ g~si þ f~si gi ð~xi Þ
ei Þ ¼ f i 1 ðXi 1 Þ f i 1 ðf i 1 ¼ f i for i ¼ 1Þ, and ðxi ~xi Þ. Since f i ð~xi Þ f i ð X r r c i i i g ðxi Þ g ðXi Þ gr , we have i
i
i
ei ; Xi Þ ¼ slððf gÞi ; X ei ; Xi Þ f i gi þ f i 1 gi : slðhi ; X s r r s ðfsi
ei ; Xi Þ f i gi þ f i gi 1 . Therefore hi ¼ Similarly, we can obtain slðhi ; X r s s r s i i 1 i i i i i 1
gr þ fr gs Þ \ ðfr gs þ fs gr Þ.
e, x 2 X, On the other hand, hi ðXi Þ ¼ f i ðXi Þ gi ðXi Þ fri gri , and for all ~x 2 X 1 e 1 e 1 e 1 we have hð~xÞ ¼ f ð~xÞ gð~xÞ 2 f ð X 1 Þ g ð X 1 Þ ¼ h ð X 1 Þ ¼ hc , hðx1 ; . . . ; xi ; ~xiþ1 ; . . . ; ~xn Þ 2 hð~xÞ þ
i X j¼1
ej ; Xj Þ ðXj X ej Þ; slðhj ; X
418
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
and ei 1 ; . . . ; X en Þ hð X eÞ þ hi ðXi Þ ¼ hðX1 ; . . . ; Xi ; X
i X
ej Þ: hjs ðXj X
j¼1
Pi
ej Þ. Thus, hir ¼ ðfri gri Þ \ ðh1c þ j¼1 hjs ðXj X For hðxÞ ¼ f ðxÞ=gðxÞ, the calculations of hic and hir can be derived similarly. ei , xi 2 Xi , there exist f~i 2 f i , g~i 2 gi with For any fixed ~xi 2 X s s s s f i ðxi Þ f i ð~xi Þ gi ðxi Þ gi ð~xi Þ ½f i ð~xi Þ þ f~si ðxi ~xi Þ gi ð~xi Þ f i ð~xi Þ ½gi ð~xi Þ þ g~si ðxi ~xi Þ gi ðxi Þ gi ð~xi Þ f~i gi ð~xi Þ f i ð~xi Þ g~si f~si hi ð~xi Þ g~si ~
ðx
ðxi xei Þ x Þ ¼ ¼ s i i gi ðxi Þ gi ð~xi Þ gi ðxi Þ
¼
ei ; Xi Þ ðf i hi gi Þ=gi . Hence, slðhi Þ ¼ slðf i =gi Þ ¼ slðhi ; X s c s r i eiþ1 ; . . . ; X en ÞÞ, i ¼ 1; . . . ; n, For h ¼ uðf Þ, since h ðyÞ ¼ uðf ðX1 ; . . . ; Xi 1 ; y; X 1 1 e 1 i i e1 ; . . . ; X en ÞÞ uðf Þ, and h ¼ h ð X ei Þ ¼ hi 1 ðXi 1 Þ we get hc ¼ h ð X 1 Þ ¼ uðf ð X c c ei , xi 2 Xi , ~xi 6¼ xi , if f i ðxi Þ 6¼ f i ð~xi Þ, we have xi 2 X hi 1 r . For any fixed ~ hi ðxi Þ hi ð~xi Þ uðf i ðxi ÞÞ uðf i ð~xi ÞÞ ¼ xi ~xi xi ~xi uðf i ðxi ÞÞ uðf i ð~xi ÞÞ f i ðxi Þ f i ð~xi Þ
¼ slðu; fci ; fri Þ fsi ; f i ðxi Þ f i ð~xi Þ xi ~xi if f i ðxi Þ ¼ f i ð~xi Þ for xi 6¼ ~xi , then 0 2 fsi , and this yields slðuðf i ÞÞ slðu; fci ; fri Þ fsi . The calculation of fri is clear. In order to get the interval extension hnr of h over X D in practical applications of Theorem 1. we must evaluate and store hi ðX i Þ, for 0 6 i 6 n 1, and only the calculation of h1c is necessary for hic . Moreover, hir is the intersection of two parts, and combines the naive interval evaluation with centered forms (see [5]). This may yield a sharper inclusion interval for the range of hi over Xi , i ¼ 1; . . . ; n. It has already been observed in [4,6] that slðu; fci ; fri Þ can be replaced by 0 u ðfci [ fri Þ when u is differentiable, where [ denotes the convex hull of the union of the intervals fci and fri . However, it is very often that slðu; fci ; fri Þ u0 ðfci [ fri Þ. In some special cases, such as for locally convex or concave functions, slðu; fci ; fri Þ can be evaluated explicitly, and may sharpen the inclusion interval for slopes significantly. The following theorem gives the form of the interval extension of the function uðf Þ, when u is a locally convex or concave function.
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
419
e D, X D. Assume that (f i ; f i ; f i ) exTheorem 2. Let f : D Rn ! R1 , X c r s ei . Let g : D b R1 ! R1 pands component functions f i of f in Xi with respect to X b for 1 6 i 6 n, where f i [ f i is the close hull of f i [ f i . Define and fci [ fri D c r c r i i hic :¼ gðfci Þ, hir :¼ gðfri Þ. Let fri :¼ ½f ir ; f r , fci :¼ ½f ic ; f c , slðg; fci ; fri Þ :¼ ½si ; si . When g is convex on fci [ fri , then (hic ; hir ; his ) expands the component function hi ei , where hi ¼ ½si ; si f i with of h ¼ gðf Þ in Xi with respect to X s s 8 i i gðf Þ gðf Þ > r c > if f ir 6¼ f ic > < f ir f ic b si ¼ g0 ðf i Þ if f ir ¼ f ic and g is differentiable on D r > > i i > gðf þeÞ gðf Þ : lim þ r r b if f ir ¼ f ic and g is non-differentiable on D e!0 e and
si ¼
8 i i gðf r Þ gðf c Þ > > > < f i f i r
0
i
i
c
g ðf r Þ > > > : lim gðf ir þeÞ gðf ir Þ e!0
i
if f r 6¼ f c
e
i i b if f r ¼ f c and g is differentiable on D i i b: if f r ¼ f c and g is non-differentiable on D
When g is concave on fci [ fri , then (hic ; hir ; his ) expands the component function hi of ei , where hi ¼ ½si ; si f i with h ¼ gðf Þ in Xi with respect to X s s 8 i i i i gðf r Þ gðf c Þ > > if f r 6¼ f c > < f ir f ic i i si ¼ g00 ðf i Þ b if f r ¼ f c and g is differentiable on D r > > i i > i i gðf r þeÞ gðf r Þ : lim b if f ¼ f and g is non-differentiable on D e!0 r
e
and
si ¼
8 i gðf Þ gðf i Þ > r c > > f i f i < r
if f ir 6¼ f ic
c
g00 ðf ir Þ > > > : lim þ gðf ir þeÞ gðf ir Þ e!0
c
e
b if f ir ¼ f ic and g is differentiable on D b: if f ir ¼ f ic and g is non-differentiable on D
Proof. When g is convex on fci [ fri , let z ¼ az1 þ ð1 aÞz2 , 0 < a < 1, z1 ; z; z2 2 fci [ fri . By convexity we get gðzÞ ¼ gðaz1 þ ð1 aÞz2 Þ 6 agðz1 Þ þ ð1 aÞgðz2 Þ and gðz1 Þ gðz2 Þ gðz1 Þ gðz2 Þ gðz1 Þ gðz2 Þ 6 ¼ ð1 aÞ z1 z2 z1 z z1 z
for z1 < z < z2 ;
hence slðg; z1 ; z2 Þ increases when z1 increases. Similarly, we have gðz1 Þ gðz2 Þ gðzÞ gðz1 Þ 6 z1 z2 z z1
for z1 < z2 < z;
420
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
which proves slðg; z1 ; z2 Þ increases when z2 increases. Thus slðg; z1 ; z2 Þ achieves its extreme value at the extreme point either fci or fri . When g is concave, the result can be similarly proved. Note that Theorem 2 yields sharper slopes for many functions, and it extends the expansion principle to non-differentiable functions. In practical computations, (i) when the expansion (fci ; fri ; fsi ) of f with e may not be exe Rn is evaluated, the interval vector X respect to some X actly representable on the computer, one possible way to solve this problem e with ~x ¼ mðX Þ; (ii) The improvement of the is to use ~x 2 Rn to replace X e Þ of the function f depends on the quality interval extension f ð~xÞ þ S ðX X of the interval vector S. According Theorems 1 and 2, the interval vector S (or the slope fs ) can be improved through intersection P of different naive ej Þ for interval evaluation (fri gri ) and centered forms h1c þ ij¼1 hjs ðXj X general functions and explicit treatment for locally convex or concave functions. To demonstrate performance using Theorems 1 and 2, as an example, we consider the function f ðx1 ; x2 Þ ¼ expðx1 x2 Þ x1 x2 for X ¼ X1 X2 ¼ e ¼X e1 X e2 ¼ fð0; 1Þg. We compare the following three ½ 1; 1 ½0; 2, X methods for expanding this function f : Method 1: gradients rf Method 2: slopes according to [4,6] Method 3: slopes as described above (using Theorems 1 and 2) e1 ; X e2 Þ, Let z1 ¼ X1 , z2 ¼ X2 , z3 ¼ z1 z2 , z4 ¼ ez3 , z5 ¼ z4 z3 , ðRzi Þ0 ¼ zi ð X > e e ðRzi Þ1 ¼ zi ðX1 ; X 2 Þ, ðRzi Þ2 ¼ zi ðX1 ; X2 Þ, Szi ¼ slðzi ; X ; X Þ ¼ ððSzi Þ1 ; ðSzi Þ2 Þ , rzi ¼ rzi ðX1 ; X2 Þ ¼ ððrzi Þ1 ; ðrzi Þ2 Þ> . Calculation results are listed in Table 1, where the final values z5 gives an interval containing the range of f on X1 X2 , i.e. z5 ðX1 ; X2 Þ is an interval extension of f on X1 X2 . By Table 1, the results of Method 2 are exactly the same for the componentwise definition of gradients (Method 1). The estimation for the range f ðX1 ; X2 Þ is the same as the estimation for Method 1. It cannot be improved by using e 1 ; X2 X e 2 Þ> g \ f r : f ðX1 ; X2 Þ ffc þ fs ðX1 X Method 3 stores more information and computes better inclusions. The last row z5 shows a sharper inclusion ½ 1:865; 7:390 for the range f ðX1 ; X2 Þ and sharper slopes ½ 0:367; 0:719 ½ 5:671; 5:671 in Method 3 as compared to Methods 1 and 2.
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
421
Table 1 Numerical results of three methods z1
z2
z3
z4
z5
Method 1
ðRzi Þ2 ðrzi Þ1 ðrzi Þ2
[)1,1] 1 0
[0,2] 0 1
[)2,2] [0,2] [)1,1]
[0.135,7.390] [0,14.779] [)7.390,7.390]
[)1.865,9.39] [)2,14.779] [)8.39,8.39]
Method 2
ðRzi Þ0 ðRzi Þ2 ðSzi Þ1 ðSzi Þ2
0 [)1,1] 1 0
1 [0,2] 0 1
0 [)2,2] 1 [)1,1]
1 [0.135,7.390] [0.135,7.390] [)7.390,7.390]
1 [)1.865,9.39] [)0.865,6.39] [)8.39,8.39]
Method 3
ðRzi Þ0 ðRzi Þ1 ðRzi Þ2 ðSzi Þ1 ðSzi Þ2
0 [)1,1] [)1,1] 1 0
1 1 [0,2] 0 1
0 [)1,1] [)2,2] 1 [)1,1]
1 [0.367,2.719] [0.135,7.390] [0.633,1.719] [)4.671,4.671]
1 [0.281,1.719] [)1.865,7.39] [)0.367,0.719] [)5.671,5.671]
3. Pruning test and global optimization algorithm In this section, we pay our attention to develop the methods for finding all solutions of global optimization problem (1) using function expansions proposed in previous section. Firstly, we give several theorems to determine regions in which it is guaranteed that there is no optimum. Secondly, form pruning tests from these theorems. Finally, give the interval pruning test algorithm for global optimization. Let f : D Rn ! R1 , X ¼ ðXi Þn1 2 IðDÞ with Xi ¼ ½ai ; bi , c ¼ ðci Þn1 2 X 0 X D. Assume that f is a known upper bound of f , i.e. f P f ¼ minx2X 0 f ðxÞ. Moreover, for 8x ¼ ðxi Þn1 2 X ; 8y 2 Xi , let f i ðyÞ :¼ f ðx1 ; . . . ; xi 1 ; y; ciþ1 ; . . . ; cn Þ with 1 6 i 6 n
ð4Þ
fsi ¼ slðf i ; ci ; Xi Þ :¼ ½fsIi ; fsSi
ð5Þ
and
dsm :¼ f f ðcÞ
n X
with 1 6 i 6 n
minðfsIj ðbj cj Þ; fsSj ðaj cj ÞÞ:
ð6Þ
j¼1;j6¼m
Theorem 3. Assume that conditions (4)–(6) hold, and there exists some m ð1 6 m 6 nÞ such that fsIm > 0. If dsm 6 0
then
pm 6 cm ; min f ðyÞ > f ; y2Y 1
ð7Þ
422
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
and if 0 6 dsm 6 ðbm cm ÞfsIm
then
c m 6 qm 6 bm ;
min f ðyÞ > f ; y2Y 2
ð8Þ
where pm ¼ cm þ dsm =fsSm , qm ¼ cm þ dsm =fsIm , Y 1 ¼ ðYi1 Þn1 , Y 2 ¼ ðYi2 Þn1 , Xi Xi if i 6¼ m if i 6¼ m 1 2 Yi ¼ Yi ¼ ðpm ; bm \ Xm if i ¼ m; ðqm ; bm if i ¼ m:
Proof. For any y 2 Y 1 ; y ¼ ðyi Þn1 . Let c ¼ ðci Þn1 2 Y 1 . Without lose generality, let Y 1 ¼ ðYj Þn1 . It follows from assumptions and definition of pm that if dsm 6 0, then pm 6 cm , and f ðyÞ 2 f ðcÞ þ
n X
slðf j ; cj ; Yj Þðyj cj Þ
j¼1 n X
f ðcÞ þ slðf m ; cm ; Ym Þðym cm Þ þ
slðf j ; cj ; Yj ÞðYj cj Þ
j¼1; j6¼m
¼ f ðcÞ þ ½fsIm ; fsSm ðym cm Þ þ
n X
½fsIj ; fsSj ðYj cj Þ:
j¼1; j6¼m
We prove result (7) for two cases. If ym P cm , since ym > pm ¼ cm þ dsm =fsS P cm þ dsm =fsI , we have f ðyÞ P f ðcÞ þ fsIm ðym cm Þ þ
n X
minðfsIj ðaj cj Þ; fsIj ðbj cj Þ; fsSj
j¼1; j6¼m
ðaj
cj Þ; fsSj
ðbj cj ÞÞ n X ¼ f ðcÞ þ fsIm ðym cm Þ þ minðfsIj ðbj cj Þ; fsSj ðaj cj ÞÞ j¼1; j6¼m
¼
fsIm
ðym cm Þ þ f
dsm
> f P f :
If ym < cm , from pm < ym < cm , we get f ðyÞ P f ðcÞ þ fsSm ðym cm Þ þ
n X
minðfsIj ðbj cj Þ; fsSj ðbj cj ÞÞ
j¼1; j6¼m
¼ fsSm ðym cm Þ þ f dsm > f P f : This completes the proof of (7). Next, we prove (8). Let c ¼ ðci Þn1 2 Y 2 . If 0 6 dsm 6 ðbm cm ÞfsIm , then we have cm 6 qm < ym . Using a similar process to the proof of (7), for any y ¼ ðyi Þn1 2 Y 2 , we have
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431 n X
f ðyÞ 2 f ðcÞ þ ½fsIm ; fsSm ðym cm Þ þ
423
½fsIj ; fsSj ðYj cj Þ:
j¼1; j6¼m
Hence f ðyÞ P f ðcÞ þ fsIm ðym cm Þ þ
n X
minðfsIj ðbj cj Þ; fsSj ðaj cj ÞÞ
j¼1; j6¼m
¼
fsIm
which gives (8).
ðym cm Þ þ f
dsm
> f P f ;
Theorem 4. Assume that conditions (4)–(6) hold, and there exists some m ð1 6 m 6 nÞ such that fsSm < 0. If 0 6 dsm 6 fsSm ðam cm Þ then
am 6 pm 6 cm ; min f ðyÞ > f ; y2Y 1
and if dsm 6 0
cm 6 qm 6 bm ; min f ðyÞ > f ;
then
y2Y 2
where pm ¼ cm þ dsm =fsSm , qm ¼ cm þ dsm =fsIm , Y 1 ¼ ðYi1 Þn1 , Y 2 ¼ ðYi2 Þn1 , Yi1 ¼
Xi ½am ; pm Þ
if i 6¼ m if i ¼ m;
Yi2 ¼
Xi ½am ; qm Þ \ Xm
if i 6¼ m if i ¼ m:
Since the proof of Theorem 4 is similar to that of Theorem 3, it is omitted here. Theorem 5. Assume that conditions (4)–(6) hold, and there exists some mð1 6 m 6 nÞ such that fsIm < 0 < fsSm . If dsm 6 0, pm 6 cm 6 qm , then miny2Y f ðyÞ > f , where pm ¼ cm þ dsm =fsSm ; qm ¼ cm þ dsm =fsIm if i 6¼ m Xi Yi ¼ ðpm ; qm Þ \ Xm if i ¼ m:
and
Y ¼ ðYi Þn1 with
Proof. From dsm 6 0 and fsIm < 0 < fsSm , we have pm 6 cm 6 qm . On the other hand, for any y ¼ ðyi Þn1 2 Y , f ðyÞ ¼ f ðcÞ þ
n X i¼1
ðf i ðyi Þ f i ðci ÞÞ ¼ f ðcÞ þ
n X
slðf i ; ci ; Yi Þðyi ci Þ:
i¼1
Since ym 2 ðpm ; qm Þ, we have the following two cases.
424
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
If ym > cm , then f ðyÞ P f ðcÞ þ fsIm ðym cm Þ þ
n X
minðfsIj ðaj cj Þ; fsIj ðbj cj Þ; fsSj
j¼1; j6¼m
þ
ðaj n X
cj Þ; fsSj
ðbj cj ÞÞ ¼ f ðcÞ þ fsIm ðym cm Þ
minðfsIj ðbj cj Þ; fsSj ðbj cj ÞÞ
j¼1; j6¼m
¼ fsIm ðym cm Þ þ f dsm > f P f : If ym 6 cm , then f ðyÞ P f ðcÞ þ fsSm ðym cm Þ þ
n X
minðfsIj ðbj cj Þ; fsSj ðbj cj ÞÞ
j¼1; j6¼m
¼
fsSm
ðym cm Þ þ f
dsm
> f P f :
Thus for any y 2 Y , miny2Y f ðyÞ > f holds and the conclusion is proved. Theorem 6. Assume that conditions (4)–(6) hold. (I) If there exists some m ð1 6 m 6 nÞ such that fsIm 6 0 < fsSm and dsm 6 0, then pm 6 cm , minw2W 1 f ðwÞ > f , where pm ¼ cm þ dsm =fsSm ; W 1 ¼ ðWi 1 Þn1 Xi if i 6¼ m Wi 1 ¼ ðpm ; bm \ Xm if i ¼ m: dsm
with
(II) Assume that there exists some m ð1 6 m 6 nÞ such that fsIm < 0 6 fsSm and 6 0, then cm 6 qm , minw2W 2 f ðwÞ > f , where qm ¼ cm þ dsm =fsIm ; W 2 ¼ ðWi 2 Þn1 Xi if i 6¼ m Wi 2 ¼ ½am ; qm Þ \ Xm if i ¼ m:
with
The proof of Theorem 6 is similar to that of Theorem 5, and is omitted. According the above theorems, we establish an interval pruning test to reject regions in which it is guaranteed that there is no optimal solutions. Let X 2 IðX 0 Þ, X ¼ ðXi Þn1 , Xi ¼ ½ai ; bi be given and choose c 2 X such that f ðcÞ > f . For i ¼ 1; 2; . . . ; n, compute dsi :¼ f f ðcÞ
n X j¼1; j6¼i
minðfsIi ðbj cj Þ; fsSj ðaj cj ÞÞ
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
425
and pi ¼ ci þ dsi =fsSi ;
qi ¼ ci þ dsi =fsIi :
Let Ti ¼ ðti ; ti Þ denote the discarded interval in Xi . Then the discarded interval Ti can be determined according the following rules. 1. If dsi 2. If dsi 3. If 4. If 5. If 6. If
fsIi > 0, then > 0. fsSi < 0, then > 0. fsIi < 0 < fsS , 0 ¼ fsIi < fsS , fsIi < fsS ¼ 0, fsIi ¼ fsS ¼ 0,
Ti ¼ ðpi ; bi for dsi 6 0, and Ti ¼ ðminðqi ; bi Þ; maxðqi ; bi ÞÞ for Ti ¼ ½ai ; qi Þ for dsi 6 0, and Ti ¼ ðminðpi ; ai Þ; maxðpi ; ai ÞÞ for then then then then
Ti Ti Ti Ti
¼ ðpi ; qi Þ for dsi < 0, and Ti ¼ ; with dsi P 0. ¼ ðpi ; bi for dsi < 0, and Ti ¼ ; with dsi P 0. ¼ ½ai ; qi Þ for dsi < 0, and Ti ¼ ; with dsi P 0. ¼ ½ai ; bi for dsi < 0, and Ti ¼ ; with dsi P 0.
Consequently discarded part of Xi is Xi \ Ti , i ¼ 1; . . . ; n. This pruning test provides the possibility to cut away a large part of the box that is currently investigated by the procedure. The rest part is left for further pruning test. Using the pruning test, at the first step the proposed global optimization algorithm discards some parts from the initial constrained box X , and generates some subboxes in which there may exist global solutions. Then the pruning test is applied onto the each the resulting subboxes to generate a list of living subboxes, and the list is denoted by Lw . Repeat the process until the living subbox is small enough to satisfy some given stopping criteria, the iteration on the subbox is terminated, and the subbox is put into the solution list Ld . Now the proposed global optimization algorithm can be described as follows. For any X ¼ ðXi Þn1 2 Lw , Xi ¼ ½ai ; bi , let AX ¼ ðai Þn1 2 Rn , BX ¼ ðbi Þn1 2 Rn , F ðX Þ ¼ ½F ðX Þ; F ðX Þ be the interval extension of f over X . f denotes the currently minimum of f found, z denotes the set currently minimum points. Interval pruning test algorithm: Step 1. Let Lw :¼ ;, Ld :¼ ; and k :¼ 0, and 0 , 1 and X 0 be given; Step 2. Let X ¼ ðXi Þn1 :¼ X 0 , c ¼ ðci Þn1 :¼ mðX 0 Þ, then calculate F ðX Þ; wðX Þ;
slðf i ; ci ; Xi Þ :¼ ½fsIi ; fsSi ;
f ¼ minðf ðAX 0 Þ; f ðBX 0 ÞÞ;
i ¼ 1; . . . ; n;
z ¼ fxjf ðxÞ ¼ f ; x ¼ AX 0 or BX 0 g;
Step 3. Let k :¼ k þ 1, c ¼ ðci Þn1 ¼ mðX Þ, 3.1 if f ðcÞ ¼ f then z ¼ fcg [ z, choose an edge of X with maximum length to bisect X along the edge to get two subboxes X 1 and X 2 : 3.2 if f ðcÞ > f then apply the interval pruning test to X , and assume that the subboxes X 1 , and/or X 2 are left;
426
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
3.3 if f ðcÞ < f then let f :¼ f ðcÞ, z :¼ fcg, update point c (for example, let c :¼ AX or c :¼ BX ) with f ðcÞ > f , and apply the interval pruning test to X , and assume that the subboxes X 1 and/or, X 2 are left; Step 4. for j :¼ 1 to 2 do if X j ¼ ; then next j; compute F ðX j Þ; wðX j Þ, and slðf i ; ci ; Xij Þ ¼ ½fsIi ; fsSi , i ¼ 1; . . . ; n; if F ðX j Þ > f , then discards the box X j and go to next j; update current minimum, let f ¼ minðf ðAX j Þ; f ðBX j Þ; f ðcÞÞ; z ¼ fxjf ðxÞ ¼ f ; x ¼ AX j ; or BX j ; or cg, if wðX j Þ < 0 or F ðX j Þ F ðX j Þ < 1 then X j is stored into Ld else X j is stored into Lw ; end for Step 5. if Lw ¼ ; then go to Step 7 else discard all X j in Lw satisfying F ðX j Þ > f , X j 2 Lw ; Step 6. if Lw ¼ ; then go to Step 7 else choose the box X j from Lw with F ðX j Þ ¼ minfF ðX m Þ; X m 2 Lw g, let X :¼ X j , then go to Step 3; Step 7. if Ld 6¼ ; then discards X j from Ld with F ðX j Þ > f ; Step 8. Output z; Ld ; k, and ½F ; f (where F ¼ minX 2Ld F ðX ÞÞ, stop. Note that in 3.2 and 3.3 of Step 3, there are three cases, (i) both X 1 and X 2 are not empty, (ii) both X 1 and X 2 are empty, (iii) one of X 1 and X 2 may be empty. 4. Numerical results To verify the computational efficiency of the proposed global optimization algorithm, a FORTRAN 90 programs is coded to implement the algorithm on several common used test functions. Numerical experiments for the following functions are implemented. Example 1. 6 Humpback camel function f ðxÞ ¼ 4x21 2:1x41 þ x61 =3 þ x1 x2 4x22 þ 4x42 ;
X ¼ fð 0:08984 . . . ; 0:71265 . . .Þg;
3 6 x1 ; x2 6 3:
f ¼ 1:031628 . . . :
Using our algorithm we obtain: global minimizers in ½ 8:984470367E 2; 8:984184265E 2 ½7:126547917E 1; 7:12656385E 1 and
½8:984184265E 2; 8:984470367E 2 ; ½ 7:12657310E 1; 7:126557084E 1
global minimum value in ½ 1:031628453528; 1:031628453481.
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
427
Example 2. 2-dimensional Shubert function (I) ( f ðxÞ ¼
5 X
)( i cos½ði þ 1Þx1 þ i
i¼1
5 X
) i cos½ði þ 1Þx2 þ i ;
i¼1
10 6 x1 ; x2 6 10: There exists 18 global minimizer points for this problem. For example, ( 7:708313 . . . ; 7:083506 . . .) and ( 0:800321 . . . ; 1:425128 . . .) and the optimal value is f ¼ 186:730908 . . . Using our algorithm we obtain 18 boxes, each of which includes one global minimizer. For example, ½ 7:708313826676376; 7:708313696080233 ; ½ 7:083506525740551; 7:083506330331995 ½ 0:8003211653744193; 0:8003210134994406 ½ 1:425128486256343; 1:425128312077582 ...... global minimum value in ½ 186:7309088310692; 186:7309088310238. Example 3. Goldstein and Price function 2
f ðxÞ ¼ ½1 þ ðx1 þ x2 þ 1Þ ð19 14x1 þ 3x21 14x2 þ 6x1 x2 þ 3x22 Þ
½30 þ ð2x1 3x2 Þ2 ð18 32x1 þ 12x21 þ 48x2 36x1 x2 þ 27x22 Þ; 2 6 x1 ; x2 6 2: X ¼ fð0; 1Þg;
f ¼ 3:0:
Using our algorithm we can obtain X and f . Example 4. 2-dimensional Shubert function (II) ( )( ) 5 5 X X j cosððj þ 1Þx1 þ jÞ j cosððj þ 1Þx2 þ jÞ f ðxÞ ¼ j¼1
j¼1 2
2
þ ðx1 þ 1:42513Þ þ ðx2 þ 0:80032Þ ;
10 6 x1 ; x2 6 10:
There exists 760 minima for this problem. The global minimum is X ¼ fð 1:42513; 0:80032Þg, f ¼ 186:7309 . . . Using our algorithm we obtain: global minimizer is ð 1:425128455269233; 8:003210604756021E 001Þ, global minimum value in ½ 186:7309088310491; 186:7309088310151.
428
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
Example 5. f ðxÞ ¼
P5
j¼1
jj cosððj þ 1Þx þ jÞj þ 5;
10 6 x 6 10.
X ¼ f5:169026 . . . ; 4:255751 . . . ; 1:114159 . . . ; 2:027433 . . . ; 7:397344 . . . ; 8:310618 . . .g;
f ¼ 6:699793 . . .
With our algorithm we obtain: global minimizers in ½5:169026041812416; 5:169026041895110;
½ 1:114159265367170; 1:114159265284479;
½8:310618695402210; 8:310618695484903; ½2:027433388206490; 2:027433388305320;
½ 4:255751918956964; 4:255751918874267; ½ 7:397344572546756; 7:397344572464068:
global minimum value in ½6:699793775630855; 6:699793776117289. Example 6. f ðxÞ ¼ jðx 1Þ=4j þ j sinðpð1 þ ðx 1Þ=4ÞÞj þ 1; 10 6 x 6 10. X ¼ f1:0g, f ¼ 1:0: With our algorithm we obtain: global minimizer in ½9:999999999999999E 01; 1:0, global minimum value in [1.0, 1.0]. Example 7 ( f ðxÞ ¼
5 X
) kj cosððk þ 1Þx1 þ kÞj þ 5
k¼1
(
5 X
) kj cosððk þ 1Þx2 þ kÞj þ 5 ;
10 6 x1 ; x2 6 10:
k¼1
This function has a large number of local minima, and among them 36 global minima, for example, ð5:169026 . . . ; 2:02743 . . .Þ and ð8:310618; 7:3973445Þ. f ¼ 44:887236 . . . Using our algorithm we can output 36 boxes, each of which includes one global minimizer. For example, ð½5:169026041762127; 5:169026041856334; ½2:027433388204543; 2:027433388257910Þ; ð½8:310618695388026; 8:310618695422186; ½ 7:397344572547674; 7:397344572510665Þ
etc., global minimum value in ½44:887236626719900; 44:887236639367390. Pn Example 8. f ðxÞ ¼ i¼1 jxi 0:5j, xi 2 ½ 5; 5, n ¼ 15. X ¼ fð0:5Þn1 g, f ¼ 0:0. With our algorithm we obtain: global minimizer is ð4:999999999563443E 001; 5:000000000291038E 001Þ151 ;
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
429
global minimum value in ½0:000000000000000E þ 0; 1:091393642127514E 10. Example 9 f ðxÞ ¼ 20 exp
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi! n 1X 0:2 jxi j exp n i¼1
xi 2 ½ 20; 30;
n 1X cosð2pxi Þ n i¼1
! þ 20;
n ¼ 4:
X ¼ fð0:0Þ41 g; f ¼ e ¼ 2:718 . . . With our algorithm we obtain: global minimizer in 0
½ 2:646977960169689E 022; B ½ 2:646977960169689E 022; B @ ½ 2:646977960169689E 022; ½ 2:646977960169689E 022;
1 3:970466940254533E 022 3:970466940254533E 022 C C; 3:970466940254533E 022 A 1:058791184067875E 021
global minimum value in ½ 2:718281828459045; 2:718281828410238. In Table 2, we give an overview of the results obtained with 0 ¼ 1 ¼ 10 10 . For each test function we list the number of algorithm iteration (NAI), the number of function evaluations (NFE), the number of interval pruning test (NPT), the number of bisections (NB), the maximal list length (MLL) ðMLL ¼ maxk jLwj j, Lwj is the list of the living subboxes at k-iteration), and the required run-time in seconds (TIME). On the other hand, the number of global minimizers for test functions is denoted by NGM, the number of output boxes which contain the global minimizers is denoted by NOB.
Table 2 Numerical results of a set of testing functions Example
TIME
NAI
MLL
NFE
NPT
NB
NOB
NGM
(1) (2) (3) (4) (5) (6) (7) (8) (9)
0.16 0.93 0.11 0.16 0.06 0.00 2.42 18.78 0.11
939 5218 2226 516 135 14 10 421 22 795 315
50 256 82 68 33 5 945 2559 286
1219 6956 2617 827 208 23 14 545 45 591 636
756 3904 2012 324 125 8 7157 0 3
183 1314 214 192 10 6 3264 22 795 315
2 18 1 1 6 1 36 1 1
2 18 1 1 6 1 36 1 1
430
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
The computational reliability of our interval pruning test algorithm can be observed from Table 2. The algorithm finds all the global optimal solutions with reasonable amount of computational effort.
5. Concluding remarks In this paper, we present an expansion scheme f ðxÞ 2 f ð~xÞ þ S ðx ~xÞ, f : Rn ! R1 for x, ~x 2 X Rn , which frequently yields sharper inclusions for S in order to compute sharper inclusions for the range of f over a domain. The emphasis of this paper is to introduce a new pruning technique based on the expansion scheme for f which can be successfully applied to interval IBBM for global optimization. Compared with the traditional interval method of global optimization, the proposed method is superior to existing similar methods using first order information in form of derivative evaluations. Under mild assumptions on global optimization problem, the main results in the paper can be extended easily to general constrained global optimization problems.
Acknowledgements The authors would like to thank Prof. Chengxian Xu of the Faculty of Science, XiÕan Jiaotong University, for many helpful remarks. This paper is supported by The National Natural Science Foundation of China (no: 69874010).
References [1] R. Horst, H. Tuy, Global optimization (deterministic approaches), Springer-Verlag, Berlin, 1996. [2] R. Horst, P.M. Pardalos, N.V. Thoai, Introduction to global optimization, Kluwer Academic Publisher, 1995. [3] H. Ratschek, J. Rokne, New computer methods for global optimization, Ellis Horwood, Chichester, 1988. [4] R. Krawczy, A. Neumaier, Interval slopes for rational functions and associated centred forms, SIAM J. Numr. Anal. 22 (1985) 604–616. [5] L.V. Kolev, Use of interval slopes for the irrational part of factorable functions, Reliable Comput. 3 (1997) 83–93. [6] Z.H. Shen, M.A. Wolfe, On interval enclosures using slope arithmetic, Appl. Math. Comput. 39 (1990) 89–105. [7] E.R. Hansen, Global optimization using interval analysis, Marcel Dekker, New York, 1992. [8] P.P. Shen, An interval algorithm for finding all global minimizers of non-smooth functions, Chinese J. Numer. Math. Appl. 21 (3) (1999) 15–20.
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
431
[9] P. Hanse, B. Jaumard, J. Xiong, Cord-slope form of TaylorÕs expansion in univariate global optimization, J. Optim. Theory Applicat. 80 (3) (1994) 441–464. [10] L.S. Zhang, M.A. Wolfe, An interval algorithm for non-differentiable global optimization, J. Appl. Math. Comput. 63 (1994) 101–122. [11] V. Nguyen Hien, S. Jean-Jacques, Computing a global optimal solution to a design centering problem, Math. Program. 53 (1992) 111–123. [12] P.P. Shen, S.Z. Yang, An interval algorithm for finding all global minimizers of a nonsmooth function of several variable, Mathematica Applicat. 14 (1) (2001) 15–21.
Applications of interval arithmetic in non-smooth global optimization Peiping Shen
a,b,*
, Kecun Zhang a, Yanjun Wang
a
a
b
Faculty of Science, Xi’an Jiaotong University, Xi’an, Shaanxi 710049, PR China Department of Mathematics, Henan Normal University, Xinxiang 453002, PR China
Abstract In this paper, an expansion scheme for general functions is presented, which can be used to achieve better enclosure of function ranges. In the context of interval branchand-bound methods for non-smooth global optimization, a pruning technique is given via the expansion. An interval pruning test is established. This pruning test offers the possibility to cut away a large part of the currently investigated box by the optimization algorithm, which can be utilized as an accelerating device similar to the monotonicity test frequently used in interval methods for smooth problems. Numerical computation shows that the proposed method is reliable and can find all solutions of global optimization problem. Ó 2002 Elsevier Inc. All rights reserved. Keywords: Global optimization; Non-smooth function; Branch-and-bound; Interval arithmetic
1. Introduction Theory and methods of global optimization can be found many applications in practice, such as the primary design of optical thin-film, the design of the electron-filter and steam power plant etc. see [11], but due to the global requirement it is not easy in comparison with theory and methods of local optimization. Horst, Tuy and Pardalos have made a progress in this field [1,2]. In these methods for global optimization, some need to compute the gradient of
*
Corresponding author. E-mail address: [email protected] (P. Shen).
0096-3003/02/$ - see front matter Ó 2002 Elsevier Inc. All rights reserved. doi:10.1016/S0096-3003(02)00417-4
414
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
function [3,4], some require only value of function [8,10,12]. They are much more welcomed by practical workers because the given algorithms only need the value of function. This paper presents a method which only compute value of function, we concentrate on the interval branch-and-bound methods (IBBM) for global optimization. IBBM is used to find the solutions of the following global optimization problems min f ðxÞ; x2X 0
f : Rn ! R1
ð1Þ
where X 0 ¼ fx 2 Rn ja0i 6 xi 6 b0i ; ði ¼ 1; . . . ; nÞg is a box and f is a continuous function over X 0 . IBBM often incorporates a monotonicity test (see [3,7]) to delete subboxes from X 0 which do not contain any global minimum point. The test uses an interval evaluation of the derivative or slope of the objective function over the current subboxes. Moreover, the performance of an IBBM is influenced by the magnitude of the interval vector S, where S satisfies expression f ðxÞ 2 f ð~xÞ þ S ðX ~xÞ for any x; ~x 2 X X 0 ;
ð2Þ
which is usually calculated as a gradient or by means of slope. In this paper, we propose rules to evaluate the expansions of non-linear functions, present a new method for computing certified enclosures of the global minimum value and all global minimum points of non-smooth objective functions subject to a box constraint by means of the expansion. The method employs a special pruning test to cut away a large part of the living subboxes, independent of whether the slope interval contains zero or not. This work has the following three advantages in comparison with existing IBBM: (i) The proposed expansion scheme can yield sharper interval vector S and the value of f over a domain than existing methods (see Table 1). (ii) The new pruning technique uses interval expansion scheme without concerning the smoothness of the objective functions. (iii) We do not need the assumption that the objective function has only finitely many local minima, which is made in [9]. This is due to the fact that the pruning test does not delete any global minimizer. Contrary, in [9] it is only proven that the function values in the deleted parts of the current interval are not smaller than the global minimum value. Moreover, the objective function in [9] is only 1-dimensional and at least three times differentiable. So the proposed method is applicable to rather general cases. The organization of the paper is as follows. In Section 2, we give expansion schemes for general functions, and the quality of enclosure for function ranges is improved through various interval intersections and special treatment for locally convex or concave functions. In Section 3, we give an interval pruning test that is used to discard regions in which there exists no global solutions and propose algorithm for finding all solutions of global optimization problem (1). Illustrative examples are presented in Section 4. The paper concludes in Section 5.
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
415
The notations of this paper are as follows: Let IðR1 Þ be the set of real intervals in R1 , for any interval Xi ¼ ½xi ; xi R1 , then IðXi Þ ¼ fYi 2 IðR1 ÞjYi Xi R1 g. The width wðXi Þ and midpoint mðXi Þ of Xi 2 IðR1 Þ are defined by wðXi Þ ¼ xi xi and mðXi Þ ¼ ðxi þ xi Þ=2. For a real interval vector (a box) X ¼ ðXi Þn1 2 IðRn Þ, where Xi ; i ¼ 1; 2; . . . ; n, are real intervals, and the width and midpoint of X are defined by wðX Þ ¼ max1 6 i 6 n wðXi Þ; mðX Þ ¼ ðmðXi ÞÞn1 2 Rn , respectively. The function F : IðRn Þ ! IðR1 Þ is called an interval extension of a function f : Rn ! R1 in X Rn , if x 2 X implies f ðxÞ 2 F ðX Þ. In other words, f ðX Þ ¼ frag ðX Þ F ðX Þ, where f ðX Þ ¼ frag ðX Þ is the range of the function f over X , that is f ðX Þ ¼ ½minx2X f ðxÞ; maxx2X f ðxÞ. It is assumed in the following that the extension functions have the isotonicity property, i.e. X Y implies F ðX Þ F ðY Þ. The global minimum value of f on X 0 is denoted by f , and the set of global minimizer points of f on X 0 by X . Our aim is to propose an algorithm to find X and f in this paper.
2. Expansion of non-linear functions In this section we will establish the rules of expansion arithmetic to yield sharper enclosure for S in (2) and sharper interval extension of f over a given domain. First, we give definitions about expansions of general non-linear functions, then two theorems which show the calculations of S and interval extension of f are presented. e ; X D be Definition 1. Let the set valued function F : D R1 ! IðR1 Þ and X 1 1 1 given. The triple ðFc ; Fr ; Fs Þ 2 IðR Þ IðR Þ IðR Þ is said to expand F in X e , if with respect to X e : 8~x 2 X
F ð~xÞ Fc ;
8x 2 X :
F ðxÞ Fr ;
e ; 8x 2 X ; 8~ 8~x 2 X y 2 F ð~xÞ; 8y 2 F ðxÞ; 9 Fes 2 Fs :
y y~ ¼ Fes ðx ~xÞ:
e is defined by Furthermore, the slope of F with respect to X and X ( ) ~ y y e Þ; x 2 X ; ~x 2 X e ; x 6¼ ~x : e ; X Þ :¼ slðF ; X y 2 F ðX Þ; y~ 2 F ð X x ~x e Þ Fc , F ðX Þ Fr , slðF ; X e ; X Þ Fs , and From Definition 1, we have F ð X e Þ; F ðX Þ; slðF ; X e ; X ÞÞ expands F in X with respect to X e: ðF ð X Definition 2. Let f : X 0 Rn ! R1 be the objective function of problem (1), and e ¼ ðX ei Þ 2 IðX 0 Þ, X ¼ ðXi Þ 2 IðX 0 Þ. For fixed X1 ; . . . ; Xi 1 ; X eiþ1 ; . . . ; X en , X n1 n1
416
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
eiþ1 ; . . . ; X en Þ is a the function f i : R1 ! IðR1 Þ with f i ðyÞ :¼ f ðX1 ; . . . ; Xi 1 ; y; X set valued function of the variable y 2 Xi and is called the component function for f . From Definition 2, we have ei Þ ¼ f ðX1 ; . . . ; Xi 1 ; X ei ; X eiþ1 ; . . . ; X en Þ; f iðX eiþ1 ; . . . ; X en Þ f i ðXi Þ ¼ f ðX1 ; . . . ; Xi 1 ; Xi ; X ei Þ; f i ðXi Þ; slðf i ; X ei ; Xi ÞÞ expands f i in Xi with respect to X ei . and ðf i ð X e For any x 2 X ; ~x 2 X , and any 0 6 i 6 n, it can be inductively obtained that i X ej ; Xj Þ ðXj X ej Þ: ð3Þ f ðx1 ; . . . ; xi ; ~xiþ1 ; ; ~xn Þ 2 f i ðXi Þ f ð~xÞ þ slðf j ; X j¼1
ei Þ ¼ f i 1 ðXi 1 Þ; i ¼ 2; . . . ; n; From Definitions 1 and 2, we have: (i) f i ð X e and (ii) for any fixed ~x 2 X , the right hand side of (3) with i ¼ n is an interval e ; X ÞÞ extension of f over X , and ðslðf i ; X n1 can be viewed as one of S in (2). In the following, we give the rule of expansion arithmetic for the evaluations of component functions in order to compute (3), for i ¼ 1; 2; . . . ; n, let e ¼ ðX ei Þ , X ¼ ðXi Þ IðRn Þ. X n1 n1 If f ðxÞ ¼ kða constant for 8x 2 X Þ, then the triple ðk; k; 0Þ expands its ei , and f i ðXi Þ ¼ k for all component function f i in Xi with respect to X i ¼ 1; . . . ; n. If f ðxÞ ¼ xj ðj ¼ 1; . . . ; nÞ for all x 2 X , then ðfci ; fri ; fsi Þ expands its compoei , where nent function f i in Xi with respect to X ( ( ej if i 6 j; ej if i < j; X X i i e i i fr ¼ f ðXi Þ ¼ and fc ¼ f ð X i Þ ¼ Xj if i > j; Xj if i P j 1 if i ¼ j; e i ; Xi Þ ¼ fsi ¼ slðf i ; X 0 if i 6¼ j: The interval extension for general functions can be constructed by recursively using the results in the following theorems. Theorems 1 and 2 show how the interval extension of a general function can be constructed from known interval extensions of simple functions. In the theorem, the interval extensions of the functions f and g are known, we give the process to construct the interval extension of the function h which is a composite function of f and g, or f and u. e D; X D. f i and gi are the comTheorem 1. Let f ; g : D Rn ! R1 and X ei Þ, f i :¼ f i ðXi Þ, ponent functions of f and g, respectively. Let fci :¼ f i ð X r ei ; Xi Þ, gi :¼ gi ð X ei Þ, gi :¼ gi ðXi Þ, gi :¼ slðgi ; X ei ; Xi Þ. If ðf i ; f i ; f i Þ and fsi ¼ slðf i ; X c r s c r s ðgci ; gri ; gsi Þ expand the component functions f i of f and gi of g in Xi with respect to ei , respectively, then ðhi ; hi ; hi Þ expands the component function hi of the function X c r s ei , where h ¼ f g for 2 fþ; ; ; =g in Xi with respect to X
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
h1c :¼ fc1 gc1 ;
hic :¼ hi 1 r
hir :¼ ðfri gri Þ \
h1c þ
i X
with i 6¼ 1
417
for 2 fþ; ; ; =g; !
ej Þ hjs ðXj X
for 2 fþ; ; ; =g;
j¼1
8 i f gi > > < s i si ðfs gr þ fri 1 gsi Þ \ ðfri gsi þ fsi gri 1 Þ; his :¼ > with gr0 :¼ gc1 and fr0 :¼ fc1 > : i ðfs hic gsi Þ=gri
for 2 fþ; g; for ¼ ; for ¼ =;
provided no division by zero occurs. If uðfcn [ frn Þ D with u : D R1 ! R1 , then for the composite function hðxÞ ¼ uðf ðxÞÞ, we have the component functions h1c :¼ uðfc1 Þ;
hic :¼ hri 1 for i 6¼ 1;
hir :¼ uðfri Þ \
h1c þ
i X
his :¼ slðu; fci ; fsi Þ fsi ; !
ej Þ : hjs ðXj X
j¼1
Proof. We only give the proof for f g, f =g, and uðf Þ, the results for other cases can be similarly proved. e1 Þ ¼ f 1 ð X e1 Þ g1 ð X e1 Þ f 1 g1 and For hðxÞ ¼ f ðxÞ gðxÞ, we have h1 ð X c c i e i 1 1 1 1 i i 1 h ð X i Þ ¼ h ðXi 1 Þ for i 6¼ 1, and hence hc ¼ fc gc , hc ¼ hr . Moreover, ( ) f i ðxi Þ gi ðxi Þ f i ð~xi Þ gi ð~xi Þ i e ei ; xi 2 Xi ; ~xi 6¼ xi : slðh ; X i ; Xi Þ ¼ ~xi 2 X xi ~xi ei ; xi 2 Xi , there exist f~i 2 f i , g~i 2 gi with For any fixed ~xi 2 X s s s s f i ðxi Þgi ðxi Þ f i ð~xi Þgi ð~xi Þ ¼ ½f i ð~xi Þ þ f~si ðxi ~xi Þ½gi ð~xi Þ þ g~si ðxi ~xi Þ f i ð~xi Þgi ð~xi Þ ¼ ½f~i gi ðxi Þ þ f i ð~xi Þ g~i ðxi ~xi Þ: s
s
Similarly, we have f ðxi Þg ðxi Þ f ð~xi Þgi ð~xi Þ ¼ ½f i ðxi Þ g~si þ f~si gi ð~xi Þ
ei Þ ¼ f i 1 ðXi 1 Þ f i 1 ðf i 1 ¼ f i for i ¼ 1Þ, and ðxi ~xi Þ. Since f i ð~xi Þ f i ð X r r c i i i g ðxi Þ g ðXi Þ gr , we have i
i
i
ei ; Xi Þ ¼ slððf gÞi ; X ei ; Xi Þ f i gi þ f i 1 gi : slðhi ; X s r r s ðfsi
ei ; Xi Þ f i gi þ f i gi 1 . Therefore hi ¼ Similarly, we can obtain slðhi ; X r s s r s i i 1 i i i i i 1
gr þ fr gs Þ \ ðfr gs þ fs gr Þ.
e, x 2 X, On the other hand, hi ðXi Þ ¼ f i ðXi Þ gi ðXi Þ fri gri , and for all ~x 2 X 1 e 1 e 1 e 1 we have hð~xÞ ¼ f ð~xÞ gð~xÞ 2 f ð X 1 Þ g ð X 1 Þ ¼ h ð X 1 Þ ¼ hc , hðx1 ; . . . ; xi ; ~xiþ1 ; . . . ; ~xn Þ 2 hð~xÞ þ
i X j¼1
ej ; Xj Þ ðXj X ej Þ; slðhj ; X
418
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
and ei 1 ; . . . ; X en Þ hð X eÞ þ hi ðXi Þ ¼ hðX1 ; . . . ; Xi ; X
i X
ej Þ: hjs ðXj X
j¼1
Pi
ej Þ. Thus, hir ¼ ðfri gri Þ \ ðh1c þ j¼1 hjs ðXj X For hðxÞ ¼ f ðxÞ=gðxÞ, the calculations of hic and hir can be derived similarly. ei , xi 2 Xi , there exist f~i 2 f i , g~i 2 gi with For any fixed ~xi 2 X s s s s f i ðxi Þ f i ð~xi Þ gi ðxi Þ gi ð~xi Þ ½f i ð~xi Þ þ f~si ðxi ~xi Þ gi ð~xi Þ f i ð~xi Þ ½gi ð~xi Þ þ g~si ðxi ~xi Þ gi ðxi Þ gi ð~xi Þ f~i gi ð~xi Þ f i ð~xi Þ g~si f~si hi ð~xi Þ g~si ~
ðx
ðxi xei Þ x Þ ¼ ¼ s i i gi ðxi Þ gi ð~xi Þ gi ðxi Þ
¼
ei ; Xi Þ ðf i hi gi Þ=gi . Hence, slðhi Þ ¼ slðf i =gi Þ ¼ slðhi ; X s c s r i eiþ1 ; . . . ; X en ÞÞ, i ¼ 1; . . . ; n, For h ¼ uðf Þ, since h ðyÞ ¼ uðf ðX1 ; . . . ; Xi 1 ; y; X 1 1 e 1 i i e1 ; . . . ; X en ÞÞ uðf Þ, and h ¼ h ð X ei Þ ¼ hi 1 ðXi 1 Þ we get hc ¼ h ð X 1 Þ ¼ uðf ð X c c ei , xi 2 Xi , ~xi 6¼ xi , if f i ðxi Þ 6¼ f i ð~xi Þ, we have xi 2 X hi 1 r . For any fixed ~ hi ðxi Þ hi ð~xi Þ uðf i ðxi ÞÞ uðf i ð~xi ÞÞ ¼ xi ~xi xi ~xi uðf i ðxi ÞÞ uðf i ð~xi ÞÞ f i ðxi Þ f i ð~xi Þ
¼ slðu; fci ; fri Þ fsi ; f i ðxi Þ f i ð~xi Þ xi ~xi if f i ðxi Þ ¼ f i ð~xi Þ for xi 6¼ ~xi , then 0 2 fsi , and this yields slðuðf i ÞÞ slðu; fci ; fri Þ fsi . The calculation of fri is clear. In order to get the interval extension hnr of h over X D in practical applications of Theorem 1. we must evaluate and store hi ðX i Þ, for 0 6 i 6 n 1, and only the calculation of h1c is necessary for hic . Moreover, hir is the intersection of two parts, and combines the naive interval evaluation with centered forms (see [5]). This may yield a sharper inclusion interval for the range of hi over Xi , i ¼ 1; . . . ; n. It has already been observed in [4,6] that slðu; fci ; fri Þ can be replaced by 0 u ðfci [ fri Þ when u is differentiable, where [ denotes the convex hull of the union of the intervals fci and fri . However, it is very often that slðu; fci ; fri Þ u0 ðfci [ fri Þ. In some special cases, such as for locally convex or concave functions, slðu; fci ; fri Þ can be evaluated explicitly, and may sharpen the inclusion interval for slopes significantly. The following theorem gives the form of the interval extension of the function uðf Þ, when u is a locally convex or concave function.
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
419
e D, X D. Assume that (f i ; f i ; f i ) exTheorem 2. Let f : D Rn ! R1 , X c r s ei . Let g : D b R1 ! R1 pands component functions f i of f in Xi with respect to X b for 1 6 i 6 n, where f i [ f i is the close hull of f i [ f i . Define and fci [ fri D c r c r i i hic :¼ gðfci Þ, hir :¼ gðfri Þ. Let fri :¼ ½f ir ; f r , fci :¼ ½f ic ; f c , slðg; fci ; fri Þ :¼ ½si ; si . When g is convex on fci [ fri , then (hic ; hir ; his ) expands the component function hi ei , where hi ¼ ½si ; si f i with of h ¼ gðf Þ in Xi with respect to X s s 8 i i gðf Þ gðf Þ > r c > if f ir 6¼ f ic > < f ir f ic b si ¼ g0 ðf i Þ if f ir ¼ f ic and g is differentiable on D r > > i i > gðf þeÞ gðf Þ : lim þ r r b if f ir ¼ f ic and g is non-differentiable on D e!0 e and
si ¼
8 i i gðf r Þ gðf c Þ > > > < f i f i r
0
i
i
c
g ðf r Þ > > > : lim gðf ir þeÞ gðf ir Þ e!0
i
if f r 6¼ f c
e
i i b if f r ¼ f c and g is differentiable on D i i b: if f r ¼ f c and g is non-differentiable on D
When g is concave on fci [ fri , then (hic ; hir ; his ) expands the component function hi of ei , where hi ¼ ½si ; si f i with h ¼ gðf Þ in Xi with respect to X s s 8 i i i i gðf r Þ gðf c Þ > > if f r 6¼ f c > < f ir f ic i i si ¼ g00 ðf i Þ b if f r ¼ f c and g is differentiable on D r > > i i > i i gðf r þeÞ gðf r Þ : lim b if f ¼ f and g is non-differentiable on D e!0 r
e
and
si ¼
8 i gðf Þ gðf i Þ > r c > > f i f i < r
if f ir 6¼ f ic
c
g00 ðf ir Þ > > > : lim þ gðf ir þeÞ gðf ir Þ e!0
c
e
b if f ir ¼ f ic and g is differentiable on D b: if f ir ¼ f ic and g is non-differentiable on D
Proof. When g is convex on fci [ fri , let z ¼ az1 þ ð1 aÞz2 , 0 < a < 1, z1 ; z; z2 2 fci [ fri . By convexity we get gðzÞ ¼ gðaz1 þ ð1 aÞz2 Þ 6 agðz1 Þ þ ð1 aÞgðz2 Þ and gðz1 Þ gðz2 Þ gðz1 Þ gðz2 Þ gðz1 Þ gðz2 Þ 6 ¼ ð1 aÞ z1 z2 z1 z z1 z
for z1 < z < z2 ;
hence slðg; z1 ; z2 Þ increases when z1 increases. Similarly, we have gðz1 Þ gðz2 Þ gðzÞ gðz1 Þ 6 z1 z2 z z1
for z1 < z2 < z;
420
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
which proves slðg; z1 ; z2 Þ increases when z2 increases. Thus slðg; z1 ; z2 Þ achieves its extreme value at the extreme point either fci or fri . When g is concave, the result can be similarly proved. Note that Theorem 2 yields sharper slopes for many functions, and it extends the expansion principle to non-differentiable functions. In practical computations, (i) when the expansion (fci ; fri ; fsi ) of f with e may not be exe Rn is evaluated, the interval vector X respect to some X actly representable on the computer, one possible way to solve this problem e with ~x ¼ mðX Þ; (ii) The improvement of the is to use ~x 2 Rn to replace X e Þ of the function f depends on the quality interval extension f ð~xÞ þ S ðX X of the interval vector S. According Theorems 1 and 2, the interval vector S (or the slope fs ) can be improved through intersection P of different naive ej Þ for interval evaluation (fri gri ) and centered forms h1c þ ij¼1 hjs ðXj X general functions and explicit treatment for locally convex or concave functions. To demonstrate performance using Theorems 1 and 2, as an example, we consider the function f ðx1 ; x2 Þ ¼ expðx1 x2 Þ x1 x2 for X ¼ X1 X2 ¼ e ¼X e1 X e2 ¼ fð0; 1Þg. We compare the following three ½ 1; 1 ½0; 2, X methods for expanding this function f : Method 1: gradients rf Method 2: slopes according to [4,6] Method 3: slopes as described above (using Theorems 1 and 2) e1 ; X e2 Þ, Let z1 ¼ X1 , z2 ¼ X2 , z3 ¼ z1 z2 , z4 ¼ ez3 , z5 ¼ z4 z3 , ðRzi Þ0 ¼ zi ð X > e e ðRzi Þ1 ¼ zi ðX1 ; X 2 Þ, ðRzi Þ2 ¼ zi ðX1 ; X2 Þ, Szi ¼ slðzi ; X ; X Þ ¼ ððSzi Þ1 ; ðSzi Þ2 Þ , rzi ¼ rzi ðX1 ; X2 Þ ¼ ððrzi Þ1 ; ðrzi Þ2 Þ> . Calculation results are listed in Table 1, where the final values z5 gives an interval containing the range of f on X1 X2 , i.e. z5 ðX1 ; X2 Þ is an interval extension of f on X1 X2 . By Table 1, the results of Method 2 are exactly the same for the componentwise definition of gradients (Method 1). The estimation for the range f ðX1 ; X2 Þ is the same as the estimation for Method 1. It cannot be improved by using e 1 ; X2 X e 2 Þ> g \ f r : f ðX1 ; X2 Þ ffc þ fs ðX1 X Method 3 stores more information and computes better inclusions. The last row z5 shows a sharper inclusion ½ 1:865; 7:390 for the range f ðX1 ; X2 Þ and sharper slopes ½ 0:367; 0:719 ½ 5:671; 5:671 in Method 3 as compared to Methods 1 and 2.
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
421
Table 1 Numerical results of three methods z1
z2
z3
z4
z5
Method 1
ðRzi Þ2 ðrzi Þ1 ðrzi Þ2
[)1,1] 1 0
[0,2] 0 1
[)2,2] [0,2] [)1,1]
[0.135,7.390] [0,14.779] [)7.390,7.390]
[)1.865,9.39] [)2,14.779] [)8.39,8.39]
Method 2
ðRzi Þ0 ðRzi Þ2 ðSzi Þ1 ðSzi Þ2
0 [)1,1] 1 0
1 [0,2] 0 1
0 [)2,2] 1 [)1,1]
1 [0.135,7.390] [0.135,7.390] [)7.390,7.390]
1 [)1.865,9.39] [)0.865,6.39] [)8.39,8.39]
Method 3
ðRzi Þ0 ðRzi Þ1 ðRzi Þ2 ðSzi Þ1 ðSzi Þ2
0 [)1,1] [)1,1] 1 0
1 1 [0,2] 0 1
0 [)1,1] [)2,2] 1 [)1,1]
1 [0.367,2.719] [0.135,7.390] [0.633,1.719] [)4.671,4.671]
1 [0.281,1.719] [)1.865,7.39] [)0.367,0.719] [)5.671,5.671]
3. Pruning test and global optimization algorithm In this section, we pay our attention to develop the methods for finding all solutions of global optimization problem (1) using function expansions proposed in previous section. Firstly, we give several theorems to determine regions in which it is guaranteed that there is no optimum. Secondly, form pruning tests from these theorems. Finally, give the interval pruning test algorithm for global optimization. Let f : D Rn ! R1 , X ¼ ðXi Þn1 2 IðDÞ with Xi ¼ ½ai ; bi , c ¼ ðci Þn1 2 X 0 X D. Assume that f is a known upper bound of f , i.e. f P f ¼ minx2X 0 f ðxÞ. Moreover, for 8x ¼ ðxi Þn1 2 X ; 8y 2 Xi , let f i ðyÞ :¼ f ðx1 ; . . . ; xi 1 ; y; ciþ1 ; . . . ; cn Þ with 1 6 i 6 n
ð4Þ
fsi ¼ slðf i ; ci ; Xi Þ :¼ ½fsIi ; fsSi
ð5Þ
and
dsm :¼ f f ðcÞ
n X
with 1 6 i 6 n
minðfsIj ðbj cj Þ; fsSj ðaj cj ÞÞ:
ð6Þ
j¼1;j6¼m
Theorem 3. Assume that conditions (4)–(6) hold, and there exists some m ð1 6 m 6 nÞ such that fsIm > 0. If dsm 6 0
then
pm 6 cm ; min f ðyÞ > f ; y2Y 1
ð7Þ
422
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
and if 0 6 dsm 6 ðbm cm ÞfsIm
then
c m 6 qm 6 bm ;
min f ðyÞ > f ; y2Y 2
ð8Þ
where pm ¼ cm þ dsm =fsSm , qm ¼ cm þ dsm =fsIm , Y 1 ¼ ðYi1 Þn1 , Y 2 ¼ ðYi2 Þn1 , Xi Xi if i 6¼ m if i 6¼ m 1 2 Yi ¼ Yi ¼ ðpm ; bm \ Xm if i ¼ m; ðqm ; bm if i ¼ m:
Proof. For any y 2 Y 1 ; y ¼ ðyi Þn1 . Let c ¼ ðci Þn1 2 Y 1 . Without lose generality, let Y 1 ¼ ðYj Þn1 . It follows from assumptions and definition of pm that if dsm 6 0, then pm 6 cm , and f ðyÞ 2 f ðcÞ þ
n X
slðf j ; cj ; Yj Þðyj cj Þ
j¼1 n X
f ðcÞ þ slðf m ; cm ; Ym Þðym cm Þ þ
slðf j ; cj ; Yj ÞðYj cj Þ
j¼1; j6¼m
¼ f ðcÞ þ ½fsIm ; fsSm ðym cm Þ þ
n X
½fsIj ; fsSj ðYj cj Þ:
j¼1; j6¼m
We prove result (7) for two cases. If ym P cm , since ym > pm ¼ cm þ dsm =fsS P cm þ dsm =fsI , we have f ðyÞ P f ðcÞ þ fsIm ðym cm Þ þ
n X
minðfsIj ðaj cj Þ; fsIj ðbj cj Þ; fsSj
j¼1; j6¼m
ðaj
cj Þ; fsSj
ðbj cj ÞÞ n X ¼ f ðcÞ þ fsIm ðym cm Þ þ minðfsIj ðbj cj Þ; fsSj ðaj cj ÞÞ j¼1; j6¼m
¼
fsIm
ðym cm Þ þ f
dsm
> f P f :
If ym < cm , from pm < ym < cm , we get f ðyÞ P f ðcÞ þ fsSm ðym cm Þ þ
n X
minðfsIj ðbj cj Þ; fsSj ðbj cj ÞÞ
j¼1; j6¼m
¼ fsSm ðym cm Þ þ f dsm > f P f : This completes the proof of (7). Next, we prove (8). Let c ¼ ðci Þn1 2 Y 2 . If 0 6 dsm 6 ðbm cm ÞfsIm , then we have cm 6 qm < ym . Using a similar process to the proof of (7), for any y ¼ ðyi Þn1 2 Y 2 , we have
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431 n X
f ðyÞ 2 f ðcÞ þ ½fsIm ; fsSm ðym cm Þ þ
423
½fsIj ; fsSj ðYj cj Þ:
j¼1; j6¼m
Hence f ðyÞ P f ðcÞ þ fsIm ðym cm Þ þ
n X
minðfsIj ðbj cj Þ; fsSj ðaj cj ÞÞ
j¼1; j6¼m
¼
fsIm
which gives (8).
ðym cm Þ þ f
dsm
> f P f ;
Theorem 4. Assume that conditions (4)–(6) hold, and there exists some m ð1 6 m 6 nÞ such that fsSm < 0. If 0 6 dsm 6 fsSm ðam cm Þ then
am 6 pm 6 cm ; min f ðyÞ > f ; y2Y 1
and if dsm 6 0
cm 6 qm 6 bm ; min f ðyÞ > f ;
then
y2Y 2
where pm ¼ cm þ dsm =fsSm , qm ¼ cm þ dsm =fsIm , Y 1 ¼ ðYi1 Þn1 , Y 2 ¼ ðYi2 Þn1 , Yi1 ¼
Xi ½am ; pm Þ
if i 6¼ m if i ¼ m;
Yi2 ¼
Xi ½am ; qm Þ \ Xm
if i 6¼ m if i ¼ m:
Since the proof of Theorem 4 is similar to that of Theorem 3, it is omitted here. Theorem 5. Assume that conditions (4)–(6) hold, and there exists some mð1 6 m 6 nÞ such that fsIm < 0 < fsSm . If dsm 6 0, pm 6 cm 6 qm , then miny2Y f ðyÞ > f , where pm ¼ cm þ dsm =fsSm ; qm ¼ cm þ dsm =fsIm if i 6¼ m Xi Yi ¼ ðpm ; qm Þ \ Xm if i ¼ m:
and
Y ¼ ðYi Þn1 with
Proof. From dsm 6 0 and fsIm < 0 < fsSm , we have pm 6 cm 6 qm . On the other hand, for any y ¼ ðyi Þn1 2 Y , f ðyÞ ¼ f ðcÞ þ
n X i¼1
ðf i ðyi Þ f i ðci ÞÞ ¼ f ðcÞ þ
n X
slðf i ; ci ; Yi Þðyi ci Þ:
i¼1
Since ym 2 ðpm ; qm Þ, we have the following two cases.
424
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
If ym > cm , then f ðyÞ P f ðcÞ þ fsIm ðym cm Þ þ
n X
minðfsIj ðaj cj Þ; fsIj ðbj cj Þ; fsSj
j¼1; j6¼m
þ
ðaj n X
cj Þ; fsSj
ðbj cj ÞÞ ¼ f ðcÞ þ fsIm ðym cm Þ
minðfsIj ðbj cj Þ; fsSj ðbj cj ÞÞ
j¼1; j6¼m
¼ fsIm ðym cm Þ þ f dsm > f P f : If ym 6 cm , then f ðyÞ P f ðcÞ þ fsSm ðym cm Þ þ
n X
minðfsIj ðbj cj Þ; fsSj ðbj cj ÞÞ
j¼1; j6¼m
¼
fsSm
ðym cm Þ þ f
dsm
> f P f :
Thus for any y 2 Y , miny2Y f ðyÞ > f holds and the conclusion is proved. Theorem 6. Assume that conditions (4)–(6) hold. (I) If there exists some m ð1 6 m 6 nÞ such that fsIm 6 0 < fsSm and dsm 6 0, then pm 6 cm , minw2W 1 f ðwÞ > f , where pm ¼ cm þ dsm =fsSm ; W 1 ¼ ðWi 1 Þn1 Xi if i 6¼ m Wi 1 ¼ ðpm ; bm \ Xm if i ¼ m: dsm
with
(II) Assume that there exists some m ð1 6 m 6 nÞ such that fsIm < 0 6 fsSm and 6 0, then cm 6 qm , minw2W 2 f ðwÞ > f , where qm ¼ cm þ dsm =fsIm ; W 2 ¼ ðWi 2 Þn1 Xi if i 6¼ m Wi 2 ¼ ½am ; qm Þ \ Xm if i ¼ m:
with
The proof of Theorem 6 is similar to that of Theorem 5, and is omitted. According the above theorems, we establish an interval pruning test to reject regions in which it is guaranteed that there is no optimal solutions. Let X 2 IðX 0 Þ, X ¼ ðXi Þn1 , Xi ¼ ½ai ; bi be given and choose c 2 X such that f ðcÞ > f . For i ¼ 1; 2; . . . ; n, compute dsi :¼ f f ðcÞ
n X j¼1; j6¼i
minðfsIi ðbj cj Þ; fsSj ðaj cj ÞÞ
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
425
and pi ¼ ci þ dsi =fsSi ;
qi ¼ ci þ dsi =fsIi :
Let Ti ¼ ðti ; ti Þ denote the discarded interval in Xi . Then the discarded interval Ti can be determined according the following rules. 1. If dsi 2. If dsi 3. If 4. If 5. If 6. If
fsIi > 0, then > 0. fsSi < 0, then > 0. fsIi < 0 < fsS , 0 ¼ fsIi < fsS , fsIi < fsS ¼ 0, fsIi ¼ fsS ¼ 0,
Ti ¼ ðpi ; bi for dsi 6 0, and Ti ¼ ðminðqi ; bi Þ; maxðqi ; bi ÞÞ for Ti ¼ ½ai ; qi Þ for dsi 6 0, and Ti ¼ ðminðpi ; ai Þ; maxðpi ; ai ÞÞ for then then then then
Ti Ti Ti Ti
¼ ðpi ; qi Þ for dsi < 0, and Ti ¼ ; with dsi P 0. ¼ ðpi ; bi for dsi < 0, and Ti ¼ ; with dsi P 0. ¼ ½ai ; qi Þ for dsi < 0, and Ti ¼ ; with dsi P 0. ¼ ½ai ; bi for dsi < 0, and Ti ¼ ; with dsi P 0.
Consequently discarded part of Xi is Xi \ Ti , i ¼ 1; . . . ; n. This pruning test provides the possibility to cut away a large part of the box that is currently investigated by the procedure. The rest part is left for further pruning test. Using the pruning test, at the first step the proposed global optimization algorithm discards some parts from the initial constrained box X , and generates some subboxes in which there may exist global solutions. Then the pruning test is applied onto the each the resulting subboxes to generate a list of living subboxes, and the list is denoted by Lw . Repeat the process until the living subbox is small enough to satisfy some given stopping criteria, the iteration on the subbox is terminated, and the subbox is put into the solution list Ld . Now the proposed global optimization algorithm can be described as follows. For any X ¼ ðXi Þn1 2 Lw , Xi ¼ ½ai ; bi , let AX ¼ ðai Þn1 2 Rn , BX ¼ ðbi Þn1 2 Rn , F ðX Þ ¼ ½F ðX Þ; F ðX Þ be the interval extension of f over X . f denotes the currently minimum of f found, z denotes the set currently minimum points. Interval pruning test algorithm: Step 1. Let Lw :¼ ;, Ld :¼ ; and k :¼ 0, and 0 , 1 and X 0 be given; Step 2. Let X ¼ ðXi Þn1 :¼ X 0 , c ¼ ðci Þn1 :¼ mðX 0 Þ, then calculate F ðX Þ; wðX Þ;
slðf i ; ci ; Xi Þ :¼ ½fsIi ; fsSi ;
f ¼ minðf ðAX 0 Þ; f ðBX 0 ÞÞ;
i ¼ 1; . . . ; n;
z ¼ fxjf ðxÞ ¼ f ; x ¼ AX 0 or BX 0 g;
Step 3. Let k :¼ k þ 1, c ¼ ðci Þn1 ¼ mðX Þ, 3.1 if f ðcÞ ¼ f then z ¼ fcg [ z, choose an edge of X with maximum length to bisect X along the edge to get two subboxes X 1 and X 2 : 3.2 if f ðcÞ > f then apply the interval pruning test to X , and assume that the subboxes X 1 , and/or X 2 are left;
426
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
3.3 if f ðcÞ < f then let f :¼ f ðcÞ, z :¼ fcg, update point c (for example, let c :¼ AX or c :¼ BX ) with f ðcÞ > f , and apply the interval pruning test to X , and assume that the subboxes X 1 and/or, X 2 are left; Step 4. for j :¼ 1 to 2 do if X j ¼ ; then next j; compute F ðX j Þ; wðX j Þ, and slðf i ; ci ; Xij Þ ¼ ½fsIi ; fsSi , i ¼ 1; . . . ; n; if F ðX j Þ > f , then discards the box X j and go to next j; update current minimum, let f ¼ minðf ðAX j Þ; f ðBX j Þ; f ðcÞÞ; z ¼ fxjf ðxÞ ¼ f ; x ¼ AX j ; or BX j ; or cg, if wðX j Þ < 0 or F ðX j Þ F ðX j Þ < 1 then X j is stored into Ld else X j is stored into Lw ; end for Step 5. if Lw ¼ ; then go to Step 7 else discard all X j in Lw satisfying F ðX j Þ > f , X j 2 Lw ; Step 6. if Lw ¼ ; then go to Step 7 else choose the box X j from Lw with F ðX j Þ ¼ minfF ðX m Þ; X m 2 Lw g, let X :¼ X j , then go to Step 3; Step 7. if Ld 6¼ ; then discards X j from Ld with F ðX j Þ > f ; Step 8. Output z; Ld ; k, and ½F ; f (where F ¼ minX 2Ld F ðX ÞÞ, stop. Note that in 3.2 and 3.3 of Step 3, there are three cases, (i) both X 1 and X 2 are not empty, (ii) both X 1 and X 2 are empty, (iii) one of X 1 and X 2 may be empty. 4. Numerical results To verify the computational efficiency of the proposed global optimization algorithm, a FORTRAN 90 programs is coded to implement the algorithm on several common used test functions. Numerical experiments for the following functions are implemented. Example 1. 6 Humpback camel function f ðxÞ ¼ 4x21 2:1x41 þ x61 =3 þ x1 x2 4x22 þ 4x42 ;
X ¼ fð 0:08984 . . . ; 0:71265 . . .Þg;
3 6 x1 ; x2 6 3:
f ¼ 1:031628 . . . :
Using our algorithm we obtain: global minimizers in ½ 8:984470367E 2; 8:984184265E 2 ½7:126547917E 1; 7:12656385E 1 and
½8:984184265E 2; 8:984470367E 2 ; ½ 7:12657310E 1; 7:126557084E 1
global minimum value in ½ 1:031628453528; 1:031628453481.
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
427
Example 2. 2-dimensional Shubert function (I) ( f ðxÞ ¼
5 X
)( i cos½ði þ 1Þx1 þ i
i¼1
5 X
) i cos½ði þ 1Þx2 þ i ;
i¼1
10 6 x1 ; x2 6 10: There exists 18 global minimizer points for this problem. For example, ( 7:708313 . . . ; 7:083506 . . .) and ( 0:800321 . . . ; 1:425128 . . .) and the optimal value is f ¼ 186:730908 . . . Using our algorithm we obtain 18 boxes, each of which includes one global minimizer. For example, ½ 7:708313826676376; 7:708313696080233 ; ½ 7:083506525740551; 7:083506330331995 ½ 0:8003211653744193; 0:8003210134994406 ½ 1:425128486256343; 1:425128312077582 ...... global minimum value in ½ 186:7309088310692; 186:7309088310238. Example 3. Goldstein and Price function 2
f ðxÞ ¼ ½1 þ ðx1 þ x2 þ 1Þ ð19 14x1 þ 3x21 14x2 þ 6x1 x2 þ 3x22 Þ
½30 þ ð2x1 3x2 Þ2 ð18 32x1 þ 12x21 þ 48x2 36x1 x2 þ 27x22 Þ; 2 6 x1 ; x2 6 2: X ¼ fð0; 1Þg;
f ¼ 3:0:
Using our algorithm we can obtain X and f . Example 4. 2-dimensional Shubert function (II) ( )( ) 5 5 X X j cosððj þ 1Þx1 þ jÞ j cosððj þ 1Þx2 þ jÞ f ðxÞ ¼ j¼1
j¼1 2
2
þ ðx1 þ 1:42513Þ þ ðx2 þ 0:80032Þ ;
10 6 x1 ; x2 6 10:
There exists 760 minima for this problem. The global minimum is X ¼ fð 1:42513; 0:80032Þg, f ¼ 186:7309 . . . Using our algorithm we obtain: global minimizer is ð 1:425128455269233; 8:003210604756021E 001Þ, global minimum value in ½ 186:7309088310491; 186:7309088310151.
428
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
Example 5. f ðxÞ ¼
P5
j¼1
jj cosððj þ 1Þx þ jÞj þ 5;
10 6 x 6 10.
X ¼ f5:169026 . . . ; 4:255751 . . . ; 1:114159 . . . ; 2:027433 . . . ; 7:397344 . . . ; 8:310618 . . .g;
f ¼ 6:699793 . . .
With our algorithm we obtain: global minimizers in ½5:169026041812416; 5:169026041895110;
½ 1:114159265367170; 1:114159265284479;
½8:310618695402210; 8:310618695484903; ½2:027433388206490; 2:027433388305320;
½ 4:255751918956964; 4:255751918874267; ½ 7:397344572546756; 7:397344572464068:
global minimum value in ½6:699793775630855; 6:699793776117289. Example 6. f ðxÞ ¼ jðx 1Þ=4j þ j sinðpð1 þ ðx 1Þ=4ÞÞj þ 1; 10 6 x 6 10. X ¼ f1:0g, f ¼ 1:0: With our algorithm we obtain: global minimizer in ½9:999999999999999E 01; 1:0, global minimum value in [1.0, 1.0]. Example 7 ( f ðxÞ ¼
5 X
) kj cosððk þ 1Þx1 þ kÞj þ 5
k¼1
(
5 X
) kj cosððk þ 1Þx2 þ kÞj þ 5 ;
10 6 x1 ; x2 6 10:
k¼1
This function has a large number of local minima, and among them 36 global minima, for example, ð5:169026 . . . ; 2:02743 . . .Þ and ð8:310618; 7:3973445Þ. f ¼ 44:887236 . . . Using our algorithm we can output 36 boxes, each of which includes one global minimizer. For example, ð½5:169026041762127; 5:169026041856334; ½2:027433388204543; 2:027433388257910Þ; ð½8:310618695388026; 8:310618695422186; ½ 7:397344572547674; 7:397344572510665Þ
etc., global minimum value in ½44:887236626719900; 44:887236639367390. Pn Example 8. f ðxÞ ¼ i¼1 jxi 0:5j, xi 2 ½ 5; 5, n ¼ 15. X ¼ fð0:5Þn1 g, f ¼ 0:0. With our algorithm we obtain: global minimizer is ð4:999999999563443E 001; 5:000000000291038E 001Þ151 ;
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
429
global minimum value in ½0:000000000000000E þ 0; 1:091393642127514E 10. Example 9 f ðxÞ ¼ 20 exp
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi! n 1X 0:2 jxi j exp n i¼1
xi 2 ½ 20; 30;
n 1X cosð2pxi Þ n i¼1
! þ 20;
n ¼ 4:
X ¼ fð0:0Þ41 g; f ¼ e ¼ 2:718 . . . With our algorithm we obtain: global minimizer in 0
½ 2:646977960169689E 022; B ½ 2:646977960169689E 022; B @ ½ 2:646977960169689E 022; ½ 2:646977960169689E 022;
1 3:970466940254533E 022 3:970466940254533E 022 C C; 3:970466940254533E 022 A 1:058791184067875E 021
global minimum value in ½ 2:718281828459045; 2:718281828410238. In Table 2, we give an overview of the results obtained with 0 ¼ 1 ¼ 10 10 . For each test function we list the number of algorithm iteration (NAI), the number of function evaluations (NFE), the number of interval pruning test (NPT), the number of bisections (NB), the maximal list length (MLL) ðMLL ¼ maxk jLwj j, Lwj is the list of the living subboxes at k-iteration), and the required run-time in seconds (TIME). On the other hand, the number of global minimizers for test functions is denoted by NGM, the number of output boxes which contain the global minimizers is denoted by NOB.
Table 2 Numerical results of a set of testing functions Example
TIME
NAI
MLL
NFE
NPT
NB
NOB
NGM
(1) (2) (3) (4) (5) (6) (7) (8) (9)
0.16 0.93 0.11 0.16 0.06 0.00 2.42 18.78 0.11
939 5218 2226 516 135 14 10 421 22 795 315
50 256 82 68 33 5 945 2559 286
1219 6956 2617 827 208 23 14 545 45 591 636
756 3904 2012 324 125 8 7157 0 3
183 1314 214 192 10 6 3264 22 795 315
2 18 1 1 6 1 36 1 1
2 18 1 1 6 1 36 1 1
430
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
The computational reliability of our interval pruning test algorithm can be observed from Table 2. The algorithm finds all the global optimal solutions with reasonable amount of computational effort.
5. Concluding remarks In this paper, we present an expansion scheme f ðxÞ 2 f ð~xÞ þ S ðx ~xÞ, f : Rn ! R1 for x, ~x 2 X Rn , which frequently yields sharper inclusions for S in order to compute sharper inclusions for the range of f over a domain. The emphasis of this paper is to introduce a new pruning technique based on the expansion scheme for f which can be successfully applied to interval IBBM for global optimization. Compared with the traditional interval method of global optimization, the proposed method is superior to existing similar methods using first order information in form of derivative evaluations. Under mild assumptions on global optimization problem, the main results in the paper can be extended easily to general constrained global optimization problems.
Acknowledgements The authors would like to thank Prof. Chengxian Xu of the Faculty of Science, XiÕan Jiaotong University, for many helpful remarks. This paper is supported by The National Natural Science Foundation of China (no: 69874010).
References [1] R. Horst, H. Tuy, Global optimization (deterministic approaches), Springer-Verlag, Berlin, 1996. [2] R. Horst, P.M. Pardalos, N.V. Thoai, Introduction to global optimization, Kluwer Academic Publisher, 1995. [3] H. Ratschek, J. Rokne, New computer methods for global optimization, Ellis Horwood, Chichester, 1988. [4] R. Krawczy, A. Neumaier, Interval slopes for rational functions and associated centred forms, SIAM J. Numr. Anal. 22 (1985) 604–616. [5] L.V. Kolev, Use of interval slopes for the irrational part of factorable functions, Reliable Comput. 3 (1997) 83–93. [6] Z.H. Shen, M.A. Wolfe, On interval enclosures using slope arithmetic, Appl. Math. Comput. 39 (1990) 89–105. [7] E.R. Hansen, Global optimization using interval analysis, Marcel Dekker, New York, 1992. [8] P.P. Shen, An interval algorithm for finding all global minimizers of non-smooth functions, Chinese J. Numer. Math. Appl. 21 (3) (1999) 15–20.
P. Shen et al. / Appl. Math. Comput. 144 (2003) 413–431
431
[9] P. Hanse, B. Jaumard, J. Xiong, Cord-slope form of TaylorÕs expansion in univariate global optimization, J. Optim. Theory Applicat. 80 (3) (1994) 441–464. [10] L.S. Zhang, M.A. Wolfe, An interval algorithm for non-differentiable global optimization, J. Appl. Math. Comput. 63 (1994) 101–122. [11] V. Nguyen Hien, S. Jean-Jacques, Computing a global optimal solution to a design centering problem, Math. Program. 53 (1992) 111–123. [12] P.P. Shen, S.Z. Yang, An interval algorithm for finding all global minimizers of a nonsmooth function of several variable, Mathematica Applicat. 14 (1) (2001) 15–21.