Research on Investigating the existence of function roots using Differential Evolution Algorithms GAO Fei, Tong Hengqing Dep. of Mathematics, Wuhan University of Technology, P.R.China, 430070 Abstract The topological degree of a function is of vital importance in investigating the existence and the number of the function’s roots which is a topic of major significance in Nonlinear Analysis. The prerequisite for the computation of the topological degree by the well known methods of Stenger and Kearfott is to resolve a sufficient refinement of the boundary of the polyhedron through Boult and Sikorski (BS) algorithm. However BS requires two linchpins, the parameter δ on the boundary of the polyhedron and an estimation of the Lipschitz constant k of the function, whose computations are a hard task in analytical way. In this paper, firstly the main concepts of differential evolution algorithm (DEA) and topological degree Theory are discussed, and then a new approach based on the DEA is proposed to calculate acceptable δ and k efficiently. The proposed approach is illustrated on several test problems and the results are reported and discussed comparing to the result existed. Key words Topological degree, Existence of Function Roots, Differential Evolution algorithm

1 Introduction Topological degree is a concept of great importance in Nonlinear Analysis. Under certain conditions, its value, computed on an open and bounded region D, is equal to the number of zeros of a function within the interior of D. Consider a class F of Lipschitz functions with constant K, defined on the unit cube C,

f : C → \ n such that for every f ∈ F we have

f ( x)

∞

≥ δ > 0 , ∀x ∈ ∂C , where ∂C is

the boundary of C. In numerous cases the parameters δ and K are not a priori known. Usually, an underestimation of δ and an overestimation of K are sufficient for the computation of the topological degree, but the computation of these estimations is a hard task and in cases where this is possible, it is performed analytically. In this paper, a new approach based on the Differential Evolution (DE) algorithm is proposed. Specifically, DE is employed to compute efficiently the parameter δ and an estimation of the Lipschitz constant K. Then, by applying Boult and Sikorski’s (BS) optimal complexity algorithm [1], a sufficient refinement of the boundary is obtained, which is required according to Stenger’s theorem [2] to compute the topological degree. DE algorithm is a novel minimization method[3][4], capable of handling non differentiable, nonlinear and multimodal objective functions, with few, easily chosen, control parameters. The crucial idea behind DE is a scheme for generating trial parameter vectors. Basically, DE adds the weighted difference between two population vectors to a third vector. This way no separate probability distribution has to be used which makes the scheme completely self-organizing[5][6]. Although the 2–dimensional case is mainly considered in this article, a generalization to higher dimensional problems is straightforward. The rest of the paper is organized as follows: in Section 2 the DE algorithm is described. In Section 3 some main concepts of the topological degree theory are briefly described and analyzed and in Section 4 experimental results are reported. The paper concludes with Section 5.

2 The Differential Evolution Algorithm The general numerical optimization problems can be defined as: Min f(x)

(1)

where x=(x1, x2,…, xn) ∈ S ⊆ R , (1≤i≤N, n ∈ ` ), xi∈[li, ui], li and ui are lower and upper values respectively, f(x) is the objective function. Suppose for a certain point x*, there exists f(x*) ≤f(x) for ∀x ∈ S , then x* and f(x*) are separately the global optimum point and its value. The solution space is defined as S0= {x|f(x)-f(x*)≤ε}, where ε is a small positive value. In order to find xi∈S0 with high probability, the typical challenges include: S0 is often very small; little a priori knowledge is available for the landscape; and calculation time is finite. Many methods based on generate-and-test have been n

proposed, such as Taboo search (TS), simulated annealing (SA), evolutionary algorithms (EAs), and others algorithms [5][6], etc. Differential Evolution (DE) algorithm grew out of Ken Price’s attempts to solve the Chebychev Polynomial fitting Problem that had been posed to him by Rainer Storn[3]. A breakthrough happened, when Ken came up with the idea of using vector differences for perturbing the vector population. The DE utilizes M n–dimensional vectors, Xi=(xi1, xi2,…, xin)∈S, i=1,…,M, as a population for each iteration, called a generation, of the algorithm. The initial population is usually taken to be uniformly distributed in the search space. At each generation, two operators, namely mutation and crossover (recombination), are applied on each individual, thus producing the new population. Then, a selection phase takes place, where each individual of the new population is compared to the corresponding individual of the old population, and the best between them is selected as a member of the population in the next generation [4]. According to the mutation operator, firstly random choose four mutually different individual in the current population X r1 , X r2 , X r3 , X r4 (G )

(G )

(G )

D

(G )

(G )

to compose a differential vector

= [ X r(1G ) − X r(2G ) ] + [ X r(3G ) − X r(4G ) ]

Then use it to add yawp to the current best individual X

(G ) best

(2)

to get a yawp vector[3][4]as below:

(G ) X βG = X best + CF × D (G )

(3)

where CF>0 is a user-defined real parameter, called mutation constant, which controls the amplification of the difference between two individuals so as to avoid search stagnation. Following the mutation phase, the crossover operator is applied on the population. For each (G )

vector X i

= ( X i(1G ) , X i(2G ) ," , X i(nG ) ) , a new vector xT( i ) = ( xT(i1) , xT(i 2) ," xT(i n) ) is generated by:

⎧⎪ X β( Gm) , if (rand (0,1) < CR) or (m = k ), = ⎨ (G ) ⎪⎩ X i m , if (rand (0,1) ≥ CR) or (m ≠ k ).

(i ) Tm

x

(4)

in the current population[4], where m=1,…,n, the index k∈{1,…,n} is randomly chosen, CR is a user-defined crossover constant[3] in the range[0,1]. In other words, the trial vector consists of some of the components of the mutant vector, and at least one of the components of a randomly selected individual of the population. Then it comes to the replacement phase. To remain the population size, we have to compare the fitness of

xT( i ) and X i(G ) , i = 1, 2,", n , then choose the better: X

( G +1) i

⎧ xT(i ) , if F ( xT(i ) ) < F ( X i( G ) ), = ⎨ (G ) otherwise. ⎩ Xi ,

(5)

Now we give the details of the DE[5] as below: Algorithm1. Differential Evolution (DE) Algorithm Step 1. Initialization. Random generate M=40 individuals in feasible space, G=0, CR=0.1, CF=0.5, Gmax=400, note the best individual in current population as

(G ) X best .

Step 2. Mutation. Random choose four mutually different individual in the current population to G

generate a yawp vector X β with (3). Step 3. Crossover. For each

X i(G ) , use (4) to get a new vector xT( i ) .

Step 4. Replacement. For each

X i( G ) and correspond xT( i ) , remain the better to the next generation

with (5). ( now )

Step 5. Updating. Find current best X best to replace the

(G ) X best if it is better than the latter, G=G +1;

Step 6. Termination condition Check. If G≤Gmax , then export the

(G ) . Else go back to Step2. X best

3 Main Concepts of Topological Degree Theory Topological degree theory has been developed as means of examining this solution set and

obtaining information on the existence of solutions, their number and their nature. Degree theory is widely used in the study of nonlinear differential (ordinary and partial) equations. It is very useful in bifurcation theory, the existence of solutions of certain partial differential equations and the use of various fixed point theorems and other aspects[7][8][9]. For the case of complex zeros the following theorem states that the total number of complex zeros can be obtained by the value of the topological degree without Picard’s extension: Theorem 1[10] Let D2∈C be an open bounded region and let f: D2→C be analytic. Suppose that f has no roots on ∂ D2 and assume that all roots of f that lie in D2 are simple. Then the total number Nr of roots of f is equal to deg[F2, D2, Θ2], where F2 ( x1 , x2 ) = ( f1 ( x1 , x2 ), f 2 ( x1 , x2 )) = (ℜ( f ( x1 + i x2 )), ℑ( f ( x1 + i x2 ))) (6)

ℜ( z ), ℑ( z ) are the real and the imaginary part of z ∈ C , respectively. Stenger’s remarkable and pioneering work[2] expresses the topological degree of a continuous mapping Fn= (f1, f2,…, fn): D n ⊂ Rn→Rn defined on a bounded domain Dn in Rn as a constant times a sum of determinants of various n×n matrices. The value of the topological degree gives information about the existence of a solution of the equation Fn(x) = Θn within Dn. The accurate computation of the topological degree of a mapping Fn at Θn relative to the bounded domain Dn using Stenger’s is heavily based on suitable assumptions, including the appropriate representation of the oriented boundary of Dn. In particular, if the boundary of Dn can be subdivided in a certain way [1][2] then Stenger’s method gives the exact value of topological degree. Otherwise, heuristic termination criteria have to be used and therefore one cannot be sure that the value of the topological degree is given correctly. Now we will concentrate on the 2–dimensional case. Definition 2[10]A segment [pi,pj] is defined to be a closed counterclockwise oriented portion of ∂ D2 with endpoints pi, and pj and interior (pi,pj). A partition P of ∂ D2 is either the empty set or a set and

{ pi }ig=1 of counterclockwise ordered points from ∂ D2 such that g

∂D2 = ∑ [ pi , pi +1 ], pg +1 = p1

(7)

i =1

Stenger[2] proved that, given a sufficient refinement[8][9]of the boundary ∂ D2 of D2, the topological degree can be computed as:

1 g deg[ F2 , D2 , Θ2 ] = ∑ (−1) ji −1 deg[ f ji +1 , [ pi , pi +1 ], 0] × sgn( f ji ( pi )) 4 i =1

(8)

where ji is the index of the component of F2=( f1, f2) that has constant sign on [pi,pi+1], f3=f1,

1 deg[ f j , [ pi , pi +1 ], 0] = sgn ⎡⎣ f j ( pi +1 ) − f j ( pi ) ⎤⎦ 2

(9)

Boult, Sikorski proposed an optimal complexity algorithm[1] for computing with certainty the topological degree for any function from a class F. For the 2–dimensional case, this class F consists of functions F2=( f1 , f2): B→R2,defined on the unit square B, satisfy the Lipschitz condition with constant K>0 and whose infinity norm along the boundary of B is at least δ>0. Then [11][12]: If K/δ≥1, then the function may have zeros in B;If the case 0.5≤K/δ≤1 it is still an open problem; If K/δ<1, then the function does not have any zeros. ∗

They established a worst-case lower bound m = 4 ⎣⎢ 4Kδ ⎦⎥ on the number of function evaluations

required to compute the topological degree of any function in the class F , by Stenger’s method[2]. If the value of the Lipschitz constant K with respect to D2 and the infinity norm δ of F2 along the boundary of D2 are known and we choose equally spaced points on the boundary of D2 separated by a

distance 1 ⎢⎣ 4Kδ ⎥⎦ , in the infinity norm, then Boult and Sikorski[1] have shown that we are able to evaluate the topological degree with certainty using Stenger’s method[2]. Thus, the values of main interest that have to be computed are δ =

min x∈∂D2

F2 ( x) ∞ , K = max x≠ y

x , y∈D2

F2 ( x) − F2 ( y ) x− y

∞

∞

.

Notice that the value of δ is always positive since the topological degree is not defined in the case where a solution of the equation F2(x)=Θ2 lies on the boundary of D2. This is also true in the Boult and Sikorski approach[1][13], since in this case the value of δ is zero and an infinite number of points have to be considered. In general, the computation of δ and K is a hard task. In the present work, DE has been employed to compute δ of F2 along the boundary of D2 as the similar effect with Particle Swarm Optimization (PSO) in[14], through subsequent minimizations on each side counterclockwise[13]Si(i=1,2,3,4) of the rectangle under consideration. Then, the smallest value obtained is considered as the value of δ. Moreover, K is estimated by repeatedly computing, through DE, the maximum of the fraction of on x, keeping y fixed, for an arbitrary large number of different y∈D2.

4 Experimental results. Comparing to traditional EA, DE algorithm has advantages in three aspects: (1) Exploration and Exploitation: In the early period of evolution, the population has great difference, so De has better capability in exploration. But as the generation goes, De has better capability in exploitation for the difference of the population is slow now. (2) Parameter Control: The DE uses crossover constant CR to control the extent of the individual parameters’ each dimension to attend the crossover phase. (3) Replacement: The replacement is done between the new and the old individuals and The DE accept the better one. This is benefit to maintain the variety of the population. So In this section, the operation of the proposed technique on different test problems is illustrated. The parameters were fixed: the size of the population was set equal to M=40, the maximum number of iterations was Gmax=400, the default values CF=0.5, CR=0.1. The individual have been constrained in the corresponding region for each test problem. The form of the test function’s quite simple, although, not differentiable, supporting the choice of DE for the minimization of the infinity norm. Example 1[13] This test problem is defined by

⎧ x 2 + x22 − 0.5 = 0, F1 = ⎨ 1 ⎩ 2 x1 x2 − 0.5 = 0,

(11)

over [0, 1] × [0,1], and it has one zero of multiplicity 2. Example 2[13] This test problem is defined by

⎧ x12 − x22 = 0, F2 = ⎨ ⎩ 2 x1 x2 = 0,

(12)

over [-1, 1] × [-1, 1], and it has two zeros. For the estimation of the Lipschitz constant K, the minimum of the infinity norm, (i.e., the value of δ), the estimation of the Lipschitz constant Kest , as well as the actual value of the Lipschitz constant, Kactual, for each test problem are reported in Table 1, where parameters resolved by PSO[14]is marked with a asterisk (*) and parameters resolved by DE. Table 1 The value of δ and K resolved by PSO [24] and DE δ Kest Kactual Test Problem δ Kest 0.5* 3.89* 0.5 3.9766 4 F1 0.8284* 3.8174* 0.8598 3.9872 4 F2 1.25 Best Fitness for Delta Average Fitness for delta

1.2

3.5

F2's Fitness for delta

F1's fitness for K

4

3 Average Fitness for K Best Fitness for K 2.5

2

1.5

1.15 1.1 1.05 1 0.95 0.9

1

0

50

100

150

200

250

300

350

400

0.85

0

50

100

150

200

250

300

350

400

Iteration Iteration Number Fig.1 The best and average Fitness for F1’s K Fig.2 The best and average Fitness for F2’s δ Fig.1 and Fig.2 shows the process the DE resolve the K of F1 and the δ of F2, where average and

best fitness denote average and best K and δ in evolution process. As the experiments show that DE is able to detect accurately the minimum of the infinity norm on the boundary of the corresponding region for each test problem, as well as to provide an acceptable estimation of the Lipschitz constant of the function. When getting acceptable δ and K as above, a sufficient refinement of the boundary can be obtained by applying Boult and Sikorski’s algorithm[1], and finally, Stenger’s theorem[2]is used to compute the topological degree with certainty.

5 Conclusions. A technique for investigating the existence of function roots has been introduced. This technique exploits the DE algorithm to compute the infinity norm on the boundary of the region, as well as to estimate the Lipschitz constant of the function under consideration. The technique has been illustrated on several test problems with satisfactory results. Further work will involve an investigation of the theoretical properties and dynamics of the proposed technique as well as applications in higher dimensional problems.

Acknowledgements F. GAO was supported by the Foundation (Grant No.XJJ2004113), Project of educational research, and the UIRT Project (Grant No. A156, A157) granted by Wuhan University of Technology inPeople's Republic of China.

References [1]

Boult T., Sikorski K., An optimal complexity algorithm for computing the topological degree in two dimensions. SIAM J. Sci. Statist. Comput., 1989,10:686-698.

[2] [3]

Stenger F., Computing the topological degree of a mapping in R . Numer. Math., 1975, 25:23-38. Storn R., Price K., Differential evolution-a simple and efficient adaptive scheme for global optimization over continuous spaces. Journal of Global Optimization, 1997, 11:341-359. Price K., An introduction to differential evolution. In D. Corne, M. Dorigo, and F. Glover, editors, New Ideas in Optimization, pages 79-108. McGraw-Hill, London (UK), 1999. Storn R. Differential evolution design of an IIR2filter. IEEE Int. Conf. on Evolutionary Computation, Nagoya, 268-273, 1996. Xie X F, Zhang W J, Bi D C. Handling equality constraints by adaptive relaxing rule for swarm algorithms. Congress on Evolutionary Computation (CEC), Oregon, USA,: 2012-2016,2004. Lloyd N.G., Degree Theory. Cambridge University Press, Cambridge, 1978. Vrahatis M.N., A short proof and a generalization of Miranda’s existence theorem. Proc. Amer. Math. Soc., 1989, 107:701–703. Vrahatis M.N., An efficient method for locating and computing periodic orbits of nonlinear mappings. J. Comp. Phys., 1995, 119:105–119. Vrahatis M.N.et al, The topological degree theory for the localization and computation of complex zeros of Bessel functions. Numer. Funct. Anal. Optimiz., 1997,18:227–234. Kearfott R.B., An efficient degree–computation method for a generalized method of bisection. Numer. Math., 32:109–127, 1979. Mourrain B., Vrahatis M.N., Yakoubsohn J.C., On the complexity of isolating real roots and computing with certainty the topological degree. J. Complexity, 2002, 18:612–640. Sikorski K.A., Optimal Solution of Nonlinear Equations. Oxford University Press, 2001. Parsopoulos, K.E., Vrahatis, M.N., Investigating the Existence of Function Roots Using Particle Swarm Optimization, Proceedings of the IEEE 2003 Congress on Evolutionary Computation, Canberra, Australia, 1448-1455, ISBN: 0-7803-7804-0, 2003.

[4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14]

n