CN106096721A - The multiple target combined schedule optimization method of Predator-prey model is thrown in based on pulse - Google Patents

Abstract

Disclosed by the invention is a kind of multiple target combined schedule optimization method throwing in Predator-prey model based on pulse, it is assumed that planted a kind of valuable medicinal plants in an ecosystem, has several with the pest population that this medicinal plants is food；There is same natural enemy populations in these pest populations, natural enemy populations constitutes Predator-prey relation with pest population；The minimizing of pest population quantity can cause the minimizing of natural enemy populations quantity；And the minimizing of natural enemy populations quantity can cause the sharp increase of pest population quantity；Periodically manually putting natural enemy populations in a suitable place to breed can make natural enemy populations quantity moment uprush, and the natural enemy populations increased suddenly can progressively suppress the sharp increase of pest population quantity；Suddenly the change of the pest population quantity that the pulse of natural enemy populations quantity increases and causes makes search can quickly jump out locally optimal solution trap；The growth change and the pulse that utilize pest population are thrown in Predator-prey model and can be quickly found out global optimum's solution scheme of multiple target combined schedule optimization problem.

Description

Based on pulse throw in predation-by food model multiple target combined schedule optimization method

Technical field

The present invention relates to intelligent optimization algorithm, be specifically related to a kind of based on pulse throw in predation-by food model multiple target Combined schedule optimization method.

Background technology

Consider that the general type of multiple target combined schedule Optimized model is as follows:

min{O₁f₁(X),O₂f₂(X),…,O_Mf_M(X)}

$\begin{array}{cc}s.t.& \left\{\begin{array}{c}{g}_{i_a}\left(X\right)\≥0,i_a\∈I\\ h_{i_b}\left(X\right)=0,i_b\∈E\\ X\∈H\⋐R^n\end{array}\right.\end{array}---\left(1\right)$

In formula:

(1)RⁿIt is that n ties up Euclidean space, the total number of variable that n is comprised by this Optimized model；

(2) X=(x₁, x₂..., x_m, x_m+1..., x_n) it is that a n ties up decision vector, wherein, front m variable x₁, x₂..., x_mIt is continuous Real-valued variable, is used for representing the resource-type parameter related in model；Rear n-m variable x_m+1..., x_nIt is 0,1 Integer type variable, is used for representing which key element in n key element must become some key elements in an optimum combination scheduling strategy, I.e. for any x_j∈{x_m+1..., x_n, if x_j=1, then it represents that jth key element is selected in the scheduling of this optimum combination Individual key element, if x_j=0, then it represents that jth key element is the most selected；Key element refers to form it is critical to of combined schedule strategy Element, such as key elements such as task, activity, manpower, goods and materials, information；

(3)f₁(X), f₂(X) ..., f_M(X) it is M object function, is used for representing M control when selecting combined schedule strategy Target call processed；

(4)O₁, O₂..., O_MFor the priority of M object function, priority order is required to meet O₁>O₂>…>O_M, i.e. mesh Scalar functions f₁(X) first require to minimize, next to that f₂(X), it is thirdly f₃(X), the like, finally require to reach That little is object function f_M(X)；

(5)Represent that key element is when selecting required satisfied i-th_aIndividual inequality constraints condition；I be inequality about The set of bundle condition number；

(6)Represent that key element is when selecting required satisfied i-th_bIndividual equality constraint；E is equality constraint bar The set of part numbering；

(7){f_i(X), i=1,2 ..., M},i_a∈I}、i_bThe mathematic(al) representation of ∈ E} does not has Restrictive condition；

(8) H is search volume, also known as solution space；

(9), when calculating, decision vector X is also referred to as trial solution；If trial solution X is unsatisfactory for constraints, then make f (X)=+ ∞。

Optimized model formula (1) be commonly used to solve the distinguished and admirable Problems of Optimal Dispatch of underground ventilating system, task scheduling optimization problem, Scheduling of resource optimization problem, job scheduling optimization problem, project scheduling optimization problem, etc..

F in Optimized model formula (1)_i(X)、Mathematic(al) representation there is no restrictive condition, traditional Mathematics Optimization Method based on continuous He the property led cannot solve this problem.

The method for solving of above-mentioned optimization problem (1) is Swarm Intelligent Algorithm, and this kind of algorithm has wide being suitable for Property.Optimized algorithm related to the present invention be competition Cooperative Evolutionary Algorithm, this algorithm by Rosin C.D., Belew R.K. at literary composition Offer " New methods for competitive coevolution, Evolutionary Computation, 1997, the Volume 5, the 1st phase, the 1-29 page " middle proposition, it is to the simulation of predation in ecogenesis: between Predator and prey Mutually producing injurious effects, either party progress of Predator and prey all can threaten the survival ability of the opposing party. The survival ability of Predator and prey is all not exclusively determined by itself, is also affected by the other side.Predator is in order to capture The survival pressure of prey can stimulate the evolution of Predator, and prey also can stimulate for the survival pressure escaping prey Prey is gradually evolved.Predator and prey mutually stimulate the other side and coevolution.

Stanley K.O., Miikkulainen R. is at document " Competitive coevolution through Evolutionary complexification, Journal ofArtificial Intelligence Research, 2004 Year, volume 21, the 63-100 page " in consider evolution complexity CompCEA algorithm carried out further patulous research；Tan T.G., Teo J., Lau H.K. is at document " Competitive coevolution with K-Random Opponents For Pareto multiobjective optimization, ICNC 2007:Third International Conference on Natural Computation, Vol 4, Proceedings, 2007:63-67 " in CompCEA is calculated Method is applied to solve Pareto multi-objective optimization question；Tan T.G., Lau H.K., Teo J. is at document " Cooperative Versus competitive coevolution for Pareto multiobjective optimization, Bio- Inspired Computational Intelligence and Applications, Springer Berlin, 2007, The 63-72 page " competition and Cooperate co-evolutionary algorithm are linked together, it is used for solving Pareto multi-objective optimization question； McIntyre A.R., Heywood M.I. is at document " Multi-objective competitive coevolution for Efficient GP classifier problem decomposition, 2007IEEE International Conference on Systems, Man and Cybernetics, Vols 1-8,2007:2582-2589 " in by multiple target Competition Cooperative Evolutionary Algorithm is used for solving GP grader resolution problem.

" multi-objective optimization question based on Ecological Prey-predator Model solves at document for Liang Wen, Luo Wenjian, Cao Xianbin, Wang Xufa Algorithm, China Science & Technology University's journal, 2005, volume 35, page 3, the 360-366 page " in according in genecology Predation model, the weights of multiple targets are corresponded to the population scale of ecosystem, utilize the dynamic week of Ecological Prey-predator Model Phase property variation characteristic, it is proposed that a kind of new Evolutionary Multiobjective Optimization.This algorithm is dynamically adjusted by Ecological Prey-predator Model The weights of each target whole, it is not necessary to artificially specify, thus decrease the dependence to problem knowledge, it is quickly found out Noninferior Solution Set.

" power plant load based on recursion multiple target Co-evolutionary Genetic Algorithm distributes at document for Wang Jianguo, Han Na, Cong Cong Research, chemical process automation and instrument, 2013, the 5th phase, the 627-630 page " in, propose to complete with Multipurpose Optimal Method Power plant load dispatching distribution, in the method, has considered the requirement of economy, rapidity and stability, has set up multiple target The mathematical model optimized, uses recursion multiple target Co-evolutionary Genetic Algorithm to solve sharing of load problem.

Guo Ning, gold paradise, Liu Wenjian are at document " manufacturing resource organization model based on virtual manufacturing cell, computer collection One-tenth manufacture system, 2010, volume 16, the 8th phase, the 649-656 page " in, it is proposed that manufacture based on virtual manufacturing cell provides Source tissue's model, has inquired into the theoretical thought of this model, gives two step construction methods of virtual manufacturing cell, and proposes profit With cooperation type multiple-objection optimization Cooperative Evolutionary Algorithm solve partner selection, construct with a plurality of processing route similarity coefficient it With the multi-objective optimization question that maximum, cost of transportation minimum and equipment produce load balancing.This multi-objective optimization question is divided Solving is multiple subproblems, evolves with many sub-populations respectively, and uses approach to cooperation between a kind of novel sub-population, finally Produce a Pareto optimal solution set.The method increase the multiformity of candidate solution, provide multiformity for production scheduling and select.

Zhang Xiaohua, Zhao Jinquan, Chen Xing Oriolus chinensis diffusus at document " containing wind energy turbine set multiple target unit compositional modeling and excellent under energy-saving and emission-reduction Change, protecting electrical power system and control, 2011, volume 39, the 17th phase, the 33-07 page " in, random for output of wind electric field Property and undulatory property, system considers upwards to be rotated down after standby and wind-powered electricity generation adds additional standby, establishes under energy-saving and emission-reduction containing wind The multiple target unit built-up pattern of electric field.By multiple targets are carried out nondimensionalization process, it is proposed that based on energy-saving and emission-reduction Decision model containing wind energy turbine set Unit Combination.This model is by regulating the weight between each target, it is possible to the energy consumption of balance system With the discharge of dusty gas, and use adaptive Cooperate co-evolutionary Algorithm for Solving.

Liu Ran, Lou Peihuang, Tang Dunbing, Yang Lei " solve the immune coevolution of U-shaped assembly line balancing scheduling problem at document Algorithm, China Mechanical Engineering, 2010, volume 21, the 7th phase, the 815-821 page ", have studied mixed flow U-shaped assembly line balancing with The multiple target integrated optimization problem of scheduling, it is proposed that a kind of multiple target clone immunity coevolution optimum based on Pareto is calculated Method.This algorithm, with two subproblems of two monoclonal antibody group's correspondence balances and scheduling, is separately encoded and coevolution, with one Polyclonal antibody group preserves optimum global solution and takes elitism strategy so that both there is cooperation between sub-population and there is also competition.Carry Go out the genotype from antibody and Phenotype evaluates antibody affinity simultaneously, and improve symbiotic partner selection mechanism to improve algorithm Constringency performance.Simulation example proves that algorithm has faster convergence rate and is more suitable for U-shaped assembling than single specie evolution algorithm Solving of line balance scheduling problem.

In sum, prior art can only solve the dimension non-combined Problems of Optimal Dispatch of the most much higher target, to dimension very High extensive solving of multiple target combined schedule optimization problem has difficulties.

Summary of the invention

In order to solve the problem that above-mentioned prior art exists, it is an object of the invention to provide a kind of input based on pulse and catch Food-and by the multiple target combined schedule optimization method of food model, i.e. MOSLO_CIPPSD method, the operator that the method constructs is permissible Fully prey between reflection different population-thrown in sky by the relation of vying each other between food relation, similar population and constant pulse Enemy population causes the unexpected change of pest population quantity, thus embodies constant pulse input Preadator prey system on multiple populations and move The basic thought of theory of mechanics；The method has global convergence.

In order to achieve the above object, the present invention adopts the following technical scheme that

A kind of based on pulse throw in predation-by food model multiple target combined schedule optimization method, i.e. MOSLO_CIPPSD Method, it is characterised in that: the general type setting multiple target combined schedule Optimized model to be solved is as follows:

min{O₁f₁(X),O₂f₂(X),…,O_Mf_M(X)}

$\begin{array}{cc}s.t.& \left\{\begin{array}{c}{g}_{i_a}\left(X\right)\≥0,i_a\∈I\\ h_{i_b}\left(X\right)=0,i_b\∈E\\ X\∈H\⋐R^n\end{array}\right.\end{array}---\left(1\right)$

In formula:

(1)RⁿIt is that n ties up Euclidean space, the total number of variable that n is comprised by this Optimized model；

(2) X=(x₁, x₂..., x_m, x_m+1..., x_n) it is that a n ties up decision vector, wherein, front m variable x₁, x₂..., x_mIt is continuous Real-valued variable, is used for representing the resource-type parameter related in model；Rear n-m variable x_m+1..., x_nIt is 0,1 Integer type variable, 0, a 1 integer type variable is also called a key element, i.e. for any x_j∈{x_m+1..., x_n, if x_j=1, Then represent jth key element selected for this optimum combination dispatch in a key element, if x_j=0, then it represents that jth key element not by Choose；

(3)f₁(X), f₂(X) ..., f_M(X) it is M object function, is used for representing M control when selecting combined schedule strategy Target call processed；

(4)O₁, O₂..., O_MFor the priority of M object function, priority order is required to meet O₁>O₂>…>O_M, i.e. mesh Scalar functions f₁(X) first require to minimize, next to that f₂(X), it is thirdly f₃(X), the like, finally require to reach That little is object function f_M(X)；

(5)Represent that key element is when selecting required satisfied i-th_aIndividual inequality constraints condition；I be inequality about The set of bundle condition number；

(6)Represent that key element is when selecting required satisfied i-th_bIndividual equality constraint；E is equality constraint bar The set of part numbering；

(7){f_i(X), i=1,2 ..., M},i_a∈I}、i_bThe mathematic(al) representation of ∈ E} does not has Restrictive condition；

(8) H is search volume, also known as solution space；

(9), when calculating, decision vector X is also referred to as trial solution；If trial solution X is unsatisfactory for constraints, then make f (X)=+ ∞；

Multiple target combined schedule Optimized model formula (1) is converted into following single goal combined schedule Optimized model:

$\begin{array}{c}\min \{F\left(X\right)=\underset{k=1}{\overset{M}{\Σ}}{O}_k{f}_k\left(X\right)\}\\ \begin{array}{cc}s.t.& \left\{\begin{array}{c}{g}_{i_a}\left(X\right)\≥0,i_a\∈I\\ h_{i_b}\left(X\right)=0,i_b\∈E\\ X\∈H\⋐R^n\end{array}\right.\end{array}\end{array}---\left(2\right)$

In formula, O_k=10^M-k；K is the numbering of object function；

The principle design of MOSLO_CIPPSD method

Assume in an ecosystem E, planted a kind of valuable medicinal plants, the species of insect pests being food with this medicinal plants Group has N kind, and they are P₁, P₂..., P_N；Each pest population has n feature, pest population P_iIt is exactly P with its character representation_i =(f_i,1,f_i,2,…,f_i,n)；In order to ensure this valuable medicinal plants natural quality, artificial dependence elements controls insect The population harm to this valuable medicinal plants；There is same natural enemy populations in this N kind pest population, this natural enemy populations is with this N kind Pest population constitutes predation-by food relation for food, i.e. natural enemy populations and pest population.

The self-seed procedure of pest population and natural enemy populations regards a continuous print process, subtracting of pest population quantity as I haven't seen you for ages causes the minimizing of natural enemy populations quantity；Owing to pest population number is a lot, the minimizing of natural enemy populations quantity can cause species of insect pests The sharp increase of group's quantity；In order to effectively control pest population, it is necessary to periodically manually put natural enemy populations in a suitable place to breed, make natural enemy populations quantity moment Uprushing, this is a kind of pulse phenomenon；Suddenly the natural enemy populations increased can progressively suppress the increasing number of pest population.

Below the solution procedure of discussion above with multiple target combined schedule Optimized model formula (2) globally optimal solution is associated Get up.

The search volume H of optimization problem is corresponding with ecosystem E, any period t, the N kind insect in this ecosystem E Population is corresponding to the N number of trial solution in the H of search volume, i.e. { X₁(t),X₂(t),…,X_N(t) }, wherein, X_i(t)=(x_i,1(t), x_i,2(t),…,x_i,n(t)), i=1,2 ..., N.Pest population P_iFeature f_i,jCorresponding to trial solution X_iIn (t) one Variable x_i,j(t)。

In summary, pest population and trial solution conceptually equivalent, the most no longer it is distinguish between.This ecosystem Each pest population in E at life span by rendering to preying on of natural enemy populations in ecosystem E, its growth conditions meeting Constantly change, this change is mapped onto the search volume H of multiple target combined schedule Optimized model formula (2), is equivalent to Trial solution transfers to another one locus from a locus；The pulse increase of natural enemy populations quantity can cause species of insect pests The unexpected change of group's quantity, this trial solution being equivalent to search volume fiercely jumps to another one position from a position, this Character is conducive to making search jump out locally optimal solution trap.

For the sake of simplicity, a locus is referred to as a state, and represents by its subscript.

Assume pest population P_iCurrent state is a, and being i.e. equivalent to location in the H of search volume is X_a.If species of insect pests Group P_iAfter being poisoned, change to new state b from current state a, be i.e. equivalent in the H of search volume from the position being presently in X_aTransfer to new position X_b.Calculate by multiple target combined schedule Optimized model formula (2), for object function F (X), if F is (X_a)>F (X_b), show new position X_bRatio original position X_aMore excellent, then it is assumed that pest population P_iEnergy for growth strong.Otherwise, if F is (X_a)≤F (X_b), show new position X_bRatio original position X_aWorse, or there is no any difference (because of new position X_bWith original position X_aTarget function value Equal, i.e. F (X_a)=F (X_b)), then it is assumed that pest population P_iEnergy for growth is weak.The pest population that energy for growth is strong, can obtain Higher probability continued growth；And the pest population that energy for growth is weak, then may stop growing.

Pest population P_iEnergy for growth power with Population Growth indices P GI (Population Growth Index, PGI) represent, pest population P_iPGI index calculation method be:

Based on above-mentioned scene, the MOSLO_ for solving multiple target combined schedule Optimized model formula (2) can be constructed CIPPSD method.

Preadator prey system kinetic model on multiple populations is thrown in constant pulse

The Preadator prey system kinetic model on multiple populations thrown in constant pulse is:

$\left\{\begin{array}{cc}\left.\begin{array}{cc}\frac{{\mathrm{dy}}_i\left(t\right)}{dt}=y_i\left(t\right)(b_i-y_i(t)-\underset{s=1,s\≠i}{\overset{N}{\Σ}}{\mathrm{\α}}_s{y}_s(t)-{\mathrm{\η}}_{i}z(t\left)\right)& i=1,2,\mathrm{...},N\\ \frac{dz\left(t\right)}{dt}=z\left(t\right)(-c+d\underset{i=1}{\overset{N}{\Σ}}{\mathrm{\η}}_i{y}_i(t\left)\right)& \end{array}\right\}& t\≠kT\\ \left.\begin{array}{cc}{\mathrm{\Δy}}_i\left(t\right)=0& i=1,2,\mathrm{...},N\\ \mathrm{\Δ}z\left(t\right)=Q& \end{array}\right\}& t=kT\end{array}\right.---\left(4\right)$

In formula: y_iT () represents t pest population P in period_iDensity, y_i(t) >=0, i=1,2 ..., N；Z (t) represents t in period The density of natural enemy populations, z (t) >=0；b_iRepresent pest population P_iRate of increase, b_i> 0, i=1,2 ..., N；η_iRepresent species of insect pests Group P_iSlip, η_i> 0, i=1,2 ..., N；α_sRepresent the contention parameter between pest population, α_s> 0, s=1,2 ..., N；d Represent the conversion ratio of natural enemy populations, d > 0；C represents the rate of increase of natural enemy populations, c > 0；T represents the input cycle of natural enemy populations；k For positive integer, k=1,2,3 ...；Q is the pulse supply volume of natural enemy populations, Q > 0.

At t in period, pest population P_iRatio shared in all pest populations is r_i(t), i=1,2 ..., N, i.e.

$r_i\left(t\right)=\frac{y_i\left(t\right)}{\underset{s=1}{\overset{N}{\Σ}}{y}_s\left(t\right)},i=1,2,...,N---\left(5\right)$

r_iT () is also called pest population P_iAccounting.Clock phase t parameter b_i, η_i, α_s, the value of c, d is respectivelyc^t, d^t；For convenience of calculating, change formula (4) into discrete recursive form, i.e.

If t ≠ kT, then

$\left\{\begin{array}{cc}{y}_i(t+1)=y_i\left(t\right)+y_i\left(t\right)(b_i^t-y_i(t)-\underset{s=1,s\≠i}{\overset{N}{\Σ}}{\mathrm{\α}}_s^t{y}_s(t)-{\mathrm{\η}}_i^{t}z(t\left)\right)& i=1,2,...,N\\ z(t+1)=z\left(t\right)+z\left(t\right)(-c^t+d^t\underset{i=1}{\overset{N}{\Σ}}{\mathrm{\η}}_i^t{y}_i(t\left)\right)& \end{array}\right.---\left(6\right)$

If t=kT, then

$\left\{\begin{array}{cc}{y}_i(t+1)=y_i\left(t\right)& i=1,2,...,N\\ z(t+1)=z\left(t\right)+Q^t& \end{array}\right.---\left(7\right)$

In formula, parameterc^t, d^tObtaining value method beb₀And b₁RepresentValue Lower limit and the upper limit, and meet b₀> 0, b₁> 0, b₀≤b₁；η₀And η₁RepresentThe lower limit of value and upper Limit, and meet η₀> 0, η₁> 0, η₀≤η₁；α₀And α₁RepresentThe lower limit of value and the upper limit, and meet α₀> 0, α₁> 0, α₀≤α₁；c^t=Rand (c₀, c₁), c₀And c₁Represent c^tThe lower limit of value and the upper limit, and meet c₀> 0, c₁> 0, c₀ ≤c₁；d^t=Rand (c₀, c₁), d₀> 0, d₀And d₁Represent d^tThe lower limit of value and the upper limit, and meet d₁> 0, d₀≤d₁；Q^t=Rand (Q₀, Q₁), Q₀And Q₁Represent Q^tThe lower limit of value and the upper limit, and meet Q₀> 0, Q₁> 0, Q₀≤Q₁；Rand (A, B) represents at [A, B] One uniform random number of interval generation, A and B is given constant, it is desirable to A≤B.

Characteristic species cluster symphysis becomes method

Period t, current pest population is P_i, characteristic species cluster symphysis becomes method as follows:

(1) producing high density pest population set AS: random choose goes out L pest population from N number of pest population, it is compiled Number formed setMake for all s ∈ { s₁, s₂..., s_L, meet r_s(t) ＞ r_i(t).L also known as For the pest population number exerted one's influence.

(2) Dominant pest population set PM is produced: first random choose goes out L population from N number of pest population, these insects The PGI index of population is than current pest population P_iPGI index high, form setWherein g₁, g₂..., g_LIt it is the numbering of these pest populations.

(3) surging pest population set SM is produced: random choose goes out L population from N number of pest population, these species of insect pestss The PGI index of group and accounting are than current population P_iPGI index and accounting high, form surging population setWherein h₁, h₂..., h_LIt it is the numbering of these pest populations；I.e. for all s ∈ { h₁, h₂..., h_L, there is PGI (X_s(t))>PGI(X_i(t)), and accounting r_s(t)>r_i(t)。

Evolutive operators

(1) bait operator.What this operator described is after pest population is preyed on by natural enemy populations, and its density declines, growth spy Levy and can change.For current pest population P_iFor, have

$v_{i,j}(t+1)=\left\{\begin{array}{cc}{r}_{i_1}\left(t\right)x_{i_1,j}\left(t\right)+r_{i_2}\left(t\right)x_{i_2,j}\left(t\right)-r_{i_3}\left(t\right)x_{i_3,j}\left(t\right)& j\≤m\\ Great(AS,i,j)& j>m\end{array}\right.---\left(8\right)$

In formula: v_i,j(t+1) it is t+1 in period current pest population P_iFeature j state value；WithRespectively t pest population in periodWithThe state value of feature j；i₁, i₂, i₃It is from { s₁, s₂..., s_LIn } with Machine is chosen out, and meets i₁≠i₂≠i₃；Great (W, j) is meant that:

k₁, k₂, k₃It is the numbering of three different populations randomly selected from set W, i.e. meets k₁≠k₂≠k₃。

(2) natural enemy operator.What this operator described is that natural enemy populations is by controlling pest population to the predation of pest population Growth, and the upgrowth situation of pest population is the change by its feature embodies.For current pest population P_iCome Say have

$v_{i,j}(t+1)=\left\{\begin{array}{cc}{r}_i\left(t\right)x_{i,j}\left(t\right)+\underset{s\∈AS}{\max }\left\{r_s\left(t\right)x_{s,j}\left(t\right)\right\}-\underset{s\∈AS}{\min }\left\{r_s\left(t\right)x_{s,j}\left(t\right)\right\}& j\≤m\\ Great(AS,i,j)& j>m\end{array}\right.---\left(9\right)$

(3) advantage operator.What this operator described is the competition between pest population, and dominant population can apply to other population Impact, i.e.

If j≤m, then

$v_{i,j}(t+1)=\left\{\begin{array}{cc}{r}_{g_a}\left(t\right)x_{g_a,j}\left(t\right)+r_{g_b}\left(t\right)x_{g_b,j}\left(t\right)-r_{g_c}\left(t\right)x_{g_c,j}\left(t\right)& \left|PM\right|\≥3\\ r_i\left(t\right)x_{i,j}\left(t\right)+r_{g_a}\left(t\right)x_{g_a,j}\left(t\right)-r_{g_b}\left(t\right)x_{g_b,j}\left(t\right)& \left|PM\right|=2\\ r_{g_a}\left(t\right)x_{g_a,j}\left(t\right)& \left|PM\right|=1\\ x_{i,j}\left(t\right)& \left|PM\right|=0\end{array}\right.---\left(10\right)$

If j > m, then

$v_{i,j}(t+1)=\left\{\begin{array}{cc}Great(PM,i,j)& \left|PM\right|\≥1\\ x_{i,j}\left(t\right)& \left|PM\right|=0\end{array}\right.---\left(11\right)$

In formula, g_a、g_b、g_cAt { g₁, g₂..., g_LRandomly choose in }, it is desirable to g_a≠g_b≠g_c。

(4) surging operator.What this operator described is the competition between pest population, and surging population can apply to other population Impact.I.e.

If j≤m, then

$v_{i,j}(t+1)=\left\{\begin{array}{cc}{r}_{h_a}\left(t\right)x_{h_a,j}\left(t\right)+r_{h_b}\left(t\right)x_{h_b,j}\left(t\right)-r_{h_c}\left(t\right)x_{h_c,j}\left(t\right)& \left|SM\right|\≥3\\ r_i\left(t\right)x_{i,j}\left(t\right)+r_{h_a}\left(t\right)x_{h_a,j}\left(t\right)-r_{h_b}\left(t\right)x_{h_b,j}\left(t\right)& \left|SM\right|=2\\ r_{h_a}\left(t\right)x_{h_a,j}\left(t\right)& \left|SM\right|=1\\ x_{i,j}\left(t\right)& \left|SM\right|=0\end{array}\right.---\left(12\right)$

If j > m, then

$v_{i,j}(t+1)=\left\{\begin{array}{cc}Great(SM,i,j)& \left|SM\right|\≥1\\ x_{i,j}\left(t\right)& \left|SM\right|=0\end{array}\right.---\left(13\right)$

In formula, h_a、h_b、h_cAt { h₁, h₂..., h_LRandomly choose in }, it is desirable to h_a≠h_b≠h_c。

(5) accretive operatos.What this operator described is the growth of population, i.e.

In formula:

X_i(t)=(x_i,1(t),x_i,2(t),…,x_i,n(t))；

V_i(t+1)=(v_i,1(t+1),v_i,2(t+1),…,v_i,n(t+1))。

The structure of MOSLO_CIPPSD method

Described MOSLO_CIPPSD method comprises the steps:

(S1) initialize:

A) t=0 is made；The all parameters related in this method are initialized by table 1；

B) initial density y of N number of pest population is determined at random₁(0), y₂(0) ..., y_N(0)；

C) initial density z (0) of natural enemy populations is determined at random；

D) N number of trial solution X is determined at random₁(0), X₂(0) ..., X_N(0)；

The obtaining value method of table 1 parameter

(S2) following operation is performed:

(S3) making t in period from 0 to G, circulation performs following step (S4)～step (S23), and wherein G is evolutionary period number；

(S4) calculate:I=1, 2 ..., N；c^t=Rand (c₀, c₁), d^t=Rand (d₀, d₁), Q^t=Rand (Q₀, Q₁)；

(S5) r is calculated by formula (5)_i(t), i=1,2 ..., N；

(S6) making i from 1 to N, circulation performs following step (S7)～step (S20)；

(S7) generate characteristic species cluster and close AS, PM, SM；

(S8) if t can not be divided exactly by T, then y is calculated by formula (6)_iAnd z (t+1) (t+1)；Otherwise, if t can be divided exactly by T, then Y is calculated by formula (7)_iAnd z (t+1) (t+1)；

(S9) making j from 1 to n, circulation performs following step (S10)～step (S18)；

(S10) calculate: p=Rand (0,1), wherein p is pest population P_iBy natural enemy populations predation and pest population P_iWith it When its similar population is vied each other, the affected actual probabilities of its growth characteristics；

(S11) if p≤E₀, then step (S12)～(S16), wherein E are performed₀For pest population because of by natural enemy populations predation and When similar population is vied each other, the affected maximum of probability of its growth characteristics；Otherwise, (S17) is gone to step；

(S12) calculate: q₀=Rand (0,1), wherein q₀For bait operator, natural enemy operator, advantage operator, surging operator quilt The actual probabilities performed；

(S13) if q₀≤ 1/4, then perform bait operator by formula (8), obtain v_i,j(t+1)；

(S14) if 1/4 < q₀≤ 1/2, then perform natural enemy operator by formula (9), obtain v_i,j(t+1)；

(S15) if 1/2 < q₀≤ 3/4, then perform advantage operator as j≤m by formula (10), obtain v_i,j(t+1)；Work as j > m Time by formula (11) perform advantage operator, obtain v_i,j(t+1)；

(S16) if 3/4 < q₀≤ 1, then as j≤m by the surging operator of formula (12), obtain v_i,j(t+1)；As j > m time by formula (13) perform surging operator, obtain v_i,j(t+1)；

(S17) if p > E₀, then v is made_i,j(t+1)=x_i,j(t)；

(S18) make j=j+1, if j≤n, then go to step (S10), otherwise go to step (S19)；

(S19) perform accretive operatos by formula (14), obtain X_i(t+1)；

(S20) make i=i+1, if i≤N, then go to step (S7), otherwise go to step (S21)；

(S21) if newly obtained globally optimal solution X^*t+1And the error between the globally optimal solution that the last time obtains meets Minimum requirements ε, then go to step (S24), otherwise goes to step (S22)；

(S22) newly obtained globally optimal solution X is preserved^*t+1；

(S23) make t=t+1, if t≤G, then turn above-mentioned steps (S4), otherwise go to step (S24)；

(S24) terminate.

Relevant parameter obtaining value method in MOSLO_PZPMD method is as shown in table 1.

Beneficial effect

Compared to the prior art the present invention, has the advantage that

1, disclosed by the invention a kind of many based on constant pulse input Preadator prey system kinetic model on multiple populations Objective cross method for optimizing scheduling, i.e. MOSLO_CIPPSD method.In the method, constant pulse is used to throw on multiple populations catching Food-by food system kinetic theory, it is assumed that in an ecosystem, plant a kind of valuable medicinal plants, with this medicinal plants Pest population for food has several；In order to ensure this valuable medicinal plants natural quality, artificial dependence elements controls The pest population harm to this valuable medicinal plants；There is same natural enemy populations in these pest populations, this natural enemy populations is with this A little pest populations constitute predation-by food relation for food, i.e. natural enemy populations and pest population.Pest population and the oneself of natural enemy populations Seed procedure is a continuous print process, and the minimizing of pest population quantity can cause the minimizing of natural enemy populations quantity；Due to insect Population number is a lot, and the minimizing of natural enemy populations quantity can cause the sharp increase of pest population quantity；In order to effectively control pest population, must Must the most manually put natural enemy populations in a suitable place to breed, make natural enemy populations quantity moment uprush, the natural enemy populations increased suddenly can progressively suppress harmful The increasing number of worm population；The growth change of pest population is equivalent to the trial solution of search volume from a position transfer to additionally One position；The pulse increase of natural enemy populations quantity can cause the unexpected change of pest population quantity, and this is equivalent to search volume Trial solution fiercely jump to another one position from a position, this character is conducive to making search to jump out locally optimal solution falling into Trap；This method has that search capability is strong and the feature of global convergence, carries for solving of multiple target combined schedule optimization problem (2) Supply a solution.

2, the search capability of MOSLO_CIPPSD method is the strongest.MOSLO_CIPPSD method includes bait operator, natural enemy Operator, advantage operator and surging operator, these operators are added significantly to its search capability.

3, model parameter value is simple.Use random method to determine based on constant pulse and throw in predation on multiple populations-by eclipse series Relevant parameter in system kinetic model and bait operator, natural enemy operator, advantage operator and the relevant parameter of surging operator, the most greatly Width decreases parameter input number, makes again model more can express practical situation.

4, the pulse increase and decrease of natural enemy populations quantity is equivalent to the trial solution of search volume and fiercely jumps to additionally from a position One position, this character is conducive to making search jump out locally optimal solution trap.

5, MOSLO_CIPPSD method take into account on multiple populations in competition process extraneous factor be discontinuously interrupted intervention Phenomenon.

6, the interaction process involved by MOSLO_CIPPSD method is rich and varied, has embodied in ecosystem common Inhomogeneity population between complicated predation-by the competitive relation between food relation and similar population.

7, the competition process involved by MOSLO_CIPPSD method has fully demonstrated rate of increase and the slip of pest population, Competitive relation between pest population, the ginseng such as cycle, pulse supply volume is thrown in the conversion ratio of natural enemy populations, rate of increase, pulse The complicated situation of change of number.

8, the feature of MOSLO_CIPPSD method of the present invention is as follows:

1) time complexity is relatively low.It is as shown in table 2 that the time complexity of MOSLO_CIPPSD method calculates process, its time Complexity and evolutionary period number G, pest population scale N, variable number n and the time complexity of each operator and other auxiliary Operation is relevant.

The time complexity computational chart of table 2MOSLO_CIPPSD method

2) MOSLO_CIPPSD method has global convergence.From bait operator, natural enemy operator, advantage operator and surging calculation The definition of son is known, generating of any one new trial solution is the most relevant with the current state of this trial solution, and with this trial solution before be The course how developing current state is unrelated, shows that the evolutionary process of MOSLO_CIPPSD method has Markov characteristic；From The definition of accretive operatos is known, the evolutionary process of MOSLO_CIPPSD method has " the poorest " characteristic；These 2 can ensure that MOSLO_CIPPSD method has global convergence, its relevant proof and document " SIS epidemic model-based Optimization, Journal of Computational Science, volume 2014,5, the 32-50 page " similar, this Bright repeat no more.

Detailed description of the invention

Below in conjunction with instantiation, the present invention is described in further detail.

(1) determine actual optimization problem to be solved, this problem is changed into multiple target combined schedule Optimized model formula (1) Described canonical form.Then, the method weighted by object function, Optimized model formula (1) is changed into single goal combination Canonical form described by Scheduling Optimization Model formula (2).

(2) parameter of MOSLO_CIPPSD method is determined by table 1.

(3) run MOSLO_CIPPSD method to solve.

(4) for following actual optimization problem, seeking n=100, overall situation when 200,400,600,800,1000,1200 is Excellent solution.

min{f₁(X),f₂(X)}

s.t.-10≤x_i≤ 10, i=1,2 ..., n-3；x_n-2+x_n-1+x_n≥1；x_n-2, x_n-1, x_n=0 or 1

$f_1\left(X\right)=\underset{i=1}{\overset{n-3}{\&Sigma;}}(x_i^2-10cos(2{\mathrm{\&pi;x}}_i)+10)+(100x_{n-2}+50x_{n-1}+x_n)$

$f_2\left(X\right)=\underset{i=1}{\overset{n-3}{\&Sigma;}}(x_i^2-20cos(2{\mathrm{\&pi;x}}_i)+20)+(150x_{n-2}+80x_{n-1}+x_n)$

A) method weighted by object function, changes into the canonical form of single-object problem by this optimization problem, I.e.

Minf (X)=10f₁(X)+f₂(X)

s.t.-10≤x_i≤ 10, i=1,2 ..., n-3；x_n-2+x_n-1+x_n≥1；x_n-2, x_n-1, x_n=0 or 1

B) parameter of MOSLO_CIPPSD method is determined by table 1, as shown in table 3.

The obtaining value method of table 3MOSLO_CIPPSD method relevant parameter

(5) using MOSLO_CIPPSD method to solve, acquired results is as shown in table 4.

Table 4 result of calculation

(6) optimal solution tried to achieve is at x_iWithin [1.126713E-8,4.231869E-8], i=1,2 ..., n-3；x_n-2 =0, x_n-1=0, x_n=1.

Claims (1)

1. based on pulse throw in predation-by food model a multiple target combined schedule optimization method, i.e. MOSLO_CIPPSD side Method, it is characterised in that: the general type setting multiple target combined schedule Optimized model to be solved is as follows:

$\begin{array}{c}\min \{O_1{f}_1\left(X\right),O_2{f}_2\left(X\right),...,O_M{f}_M\left(X\right)\}\\ \begin{array}{cc}s.t.& \left\{\begin{array}{c}{g}_{i_a}\left(X\right)\&GreaterEqual;0,i_a\&Element;I\\ h_{i_b}\left(X\right)=0,i_b\&Element;E\\ X\&Element;H\&Subset;R^n\end{array}\right.\end{array}\end{array}---\left(1\right)$

In formula:

(1)RⁿIt is that n ties up Euclidean space, the total number of variable that n is comprised by this Optimized model；

(2) X=(x₁, x₂..., x_m, x_m+1..., x_n) it is that a n ties up decision vector, wherein, front m variable x₁, x₂..., x_mIt is even Continuous Real-valued variable, is used for representing the resource-type parameter related in model；Rear n-m variable x_m+1..., x_nIt it is 0,1 integer type Variable, 0, a 1 integer type variable is also called a key element, i.e. for any x_j∈{x_m+1..., x_n, if x_j=1, then it represents that Jth key element selected for this optimum combination dispatch in a key element, if x_j=0, then it represents that jth key element is the most selected；

(3)f₁(X), f₂(X) ..., f_M(X) it is M object function, is used for M the control mesh represented when selecting combined schedule strategy Mark requirement；

(4)O₁, O₂..., O_MFor the priority of M object function, priority order is required to meet O₁>O₂>…>O_M, i.e. target letter Number f₁(X) first require to minimize, next to that f₂(X), it is thirdly f₃(X), the like, finally requirement minimizes It is object function f_M(X)；

(5)Represent that key element is when selecting required satisfied i-th_aIndividual inequality constraints condition；I is inequality constraints bar The set of part numbering；

(6)Represent that key element is when selecting required satisfied i-th_bIndividual equality constraint；E is that equality constraint is compiled Number set；

(7){f_i(X), i=1,2 ..., M},Mathematic(al) representation do not have Restrictive condition；

(8) H is search volume, also known as solution space；

(9), when calculating, decision vector X is also referred to as trial solution；If trial solution X is unsatisfactory for constraints, then make f (X)=+ ∞；

Multiple target combined schedule Optimized model formula (1) is converted into following single goal combined schedule Optimized model:

$\begin{array}{c}\min \{F\left(X\right)=\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}{O}_k{f}_k\left(X\right)\}\\ \begin{array}{cc}s.t.& \left\{\begin{array}{c}{g}_{i_a}\left(X\right)\&GreaterEqual;0,i_a\&Element;I\\ h_{i_b}\left(X\right)=0,i_b\&Element;E\\ X\&Element;H\&Subset;R^n\end{array}\right.\end{array}\end{array}---\left(2\right)$

In formula, O_k=10^M-k；K is the numbering of object function；

Described MOSLO_CIPPSD method uses constant pulse to throw in Preadator prey system kinetic theory on multiple populations, it is assumed that One ecosystem is planted a kind of valuable medicinal plants, has had several with the pest population that this medicinal plants is food；In order to Guaranteeing this valuable medicinal plants natural quality, artificial dependence elements controls the pest population danger to this valuable medicinal plants Evil；There is same natural enemy populations in these pest populations, this natural enemy populations is with these pest populations for food, i.e. natural enemy populations and evil The predation of worm population composing-by food relation；The self-seed procedure of pest population and natural enemy populations is a continuous print process, insect The minimizing of population quantity can cause the minimizing of natural enemy populations quantity；Owing to pest population number is a lot, the minimizing of natural enemy populations quantity The sharp increase of pest population quantity can be caused；In order to effectively control pest population, it is necessary to periodically manually put natural enemy populations in a suitable place to breed, make natural enemy Population quantity moment uprushes, and the natural enemy populations increased suddenly can progressively suppress the increasing number of pest population；The life of pest population Long change is equivalent to the trial solution of search volume from a position transfer to another one position；The pulse of natural enemy populations quantity increases Adding the unexpected change that can cause pest population quantity, this trial solution being equivalent to search volume fiercely jumps to additionally from a position One position, this character is conducive to making search jump out locally optimal solution trap；

Period t, pest population P_iEnergy for growth power Population Growth indices P GI represent, pest population P_iPGI index Computational methods are:

In formula, X_iT () is t pest population P in period_iCorresponding trial solution；N is insect population number；I represents pest population P_i's Numbering；

Described MOSLO_CIPPSD method comprises the steps:

(S1) initialize:

A) t=0 in period is made；The all parameters related in this method are initialized by table 1；

B) initial density y of N number of pest population is determined at random₁(0), y₂(0) ..., y_N(0)；

C) initial density z (0) of natural enemy populations is determined at random；

D) N number of trial solution X is determined at random₁(0), X₂(0) ..., X_N(0)；

The obtaining value method of table 1 parameter

(S2) following operation is performed:

(S3) making t in period from 0 to G, circulation performs following step (S4)～step (S23), and wherein G is evolutionary period number；

(S4) calculate: c^t=Rand (c₀, c₁), d^t=Rand (d₀, d₁), Q^t=Rand (Q₀, Q₁)；

In formula:c^t, d^t, Q^tIt is respectively parameter b_i, η_i, α_i, c, d, Q are in the value of t in period；b_iRepresent species of insect pests Group P_iRate of increase, b_i> 0, b₀And b₁RepresentThe lower limit of value and the upper limit, and meet b₀> 0, b₁> 0, b₀≤b₁；η_iRepresent insect Population P_iSlip, η_i> 0, η₀And η₁RepresentThe lower limit of value and the upper limit, and meet η₀> 0, η₁> 0, η₀≤η₁；α_iRepresent evil Contention parameter between worm population, α_i> 0, α₀And α₁RepresentThe lower limit of value and the upper limit, and meet α₀> 0, α₁> 0, α₀≤α₁；d Represent the conversion ratio of natural enemy populations, d > 0, d₀And d₁Represent d^tThe lower limit of value and the upper limit, and meet d₀> 0, d₁> 0, d₀≤d₁；c Represent the rate of increase of natural enemy populations, c > 0, c₀And c₁Represent c^tThe lower limit of value and the upper limit, and meet c₀> 0, c₁> 0, c₀≤c₁；Q Represent the injected volume of natural enemy populations, Q > 0, Q₀And Q₁Represent Q^tThe lower limit of value and the upper limit, and meet Q₀> 0, Q₁> 0, Q₀≤Q₁； Rand (A, B) represents that A and B is given constant, it is desirable to A≤B at one uniform random number of [A, B] interval generation；

(S5) r is calculated by formula (5)_i(t):

$r_i\left(t\right)=\frac{y_i\left(t\right)}{\underset{s=1}{\overset{N}{\mathrm{\&Sigma;}}}{y}_s\left(t\right)},i=1,2,...,N---\left(5\right)$

In formula, r_iT () is t pest population P in period_iRatio shared in all pest populations, r_iT () is also called pest population P_i Accounting；y_i(t), y_sT () represents t pest population P in period respectively_iAnd P_sDensity, y_i(t) >=0, y_s(t)≥0；

(S6) making i from 1 to N, circulation performs following step (S7)～step (S20)；

(S7) generate characteristic species cluster and close AS, PM, SM；Period, t, was P for current pest population_i, characteristic species cluster close AS, PM, The generation method of SM is as follows:

A) high density pest population set AS is produced: random choose goes out L pest population from N number of pest population, its numbering shape Become setMake for all s ∈ { s₁, s₂..., s_L, meet r_s(t) ＞ r_i(t)；L is also called and executes Add the pest population number of impact；

B) Dominant pest population set PM is produced: first random choose goes out L population from N number of pest population, these pest populations PGI index than current pest population P_iPGI index high, form setWherein g₁, g₂..., g_LIt it is the numbering of these pest populations；

C) surging pest population set SM is produced: random choose goes out L population from N number of pest population, these pest populations PGI index and accounting are than current population P_iPGI index and accounting high, form surging population setWherein h₁, h₂..., h_LIt it is the numbering of these pest populations；I.e. for all s ∈ { h₁, h₂..., h_L, there is PGI (X_s(t))>PGI(X_i(t)), and accounting r_s(t)>r_i(t)；

(S8) if t can not be divided exactly by T, then pest population P is calculated by formula (6)_iDensity y_iAnd the density z (t of natural enemy populations (t+1) +1)；Otherwise, if t can be divided exactly by T, then calculate y by formula (7)_iAnd z (t+1) (t+1)；

$\left\{\begin{array}{cc}{y}_i(t+1)=y_i\left(t\right)+y_i\left(t\right)(b_i^t-y_i(t)-\underset{s=1,s\&NotEqual;i}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_s^t{y}_s(t)-{\mathrm{\&eta;}}_i^{t}z(t\left)\right)& i=1,2,...,N\\ z(t+1)=z\left(t\right)+z\left(t\right)(-c^t+d^t\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&eta;}}_i^t{y}_i(t\left)\right)& \end{array}\right.---\left(6\right)$

$\left\{\begin{array}{cc}{y}_i(t+1)=y_i\left(t\right)& i=1,2,...,N\\ z(t+1)=z\left(t\right)+Q^t& \end{array}\right.---\left(7\right)$

In formula, T represents the input cycle of natural enemy populations；Z (t) represents the density of t natural enemy populations in period, z (t) >=0；

Formula (6), formula (7) come from the Preadator prey system kinetics on multiple populations thrown in constant pulse described by formula (4) Model:

$\left\{\begin{array}{cc}\left.\begin{array}{cc}\frac{{\mathrm{dy}}_i\left(t\right)}{dt}=y_i\left(t\right)(b_i-y_i(t)-\underset{s=1,s\&NotEqual;i}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_s{y}_s(t)-{\mathrm{\&eta;}}_{i}z(t\left)\right)& i=1,2,...,N\\ \frac{dz\left(t\right)}{dt}=z\left(t\right)(-c+d\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&eta;}}_i{y}_i(t\left)\right)& \end{array}\right\}& t\&NotEqual;kT\\ \left.\begin{array}{cc}{\mathrm{\&Delta;y}}_i\left(t\right)=0& i=1,2,...,N\\ \mathrm{\&Delta;}z\left(t\right)=Q& \end{array}\right\}& t=kT\end{array}\right.---\left(4\right)$

(S9) making j from 1 to n, circulation performs following step (S10)～step (S18)；

(S10) calculate: p=Rand (0,1), wherein p is pest population P_iBy natural enemy populations predation and pest population P_iSame with other When class population is vied each other, the affected actual probabilities of its growth characteristics；

(S11) if p≤E₀, then step (S12)～(S16), wherein E are performed₀For pest population because being preyed on and similar by natural enemy populations When population is vied each other, the affected maximum of probability of its growth characteristics；Otherwise, (S17) is gone to step；

(S12) calculate: q₀=Rand (0,1), wherein q₀It is performed for bait operator, natural enemy operator, advantage operator, surging operator Actual probabilities；

(S13) if q₀≤ 1/4, then perform bait operator by formula (8), obtain v_i,j(t+1)；

$v_{i,j}(t+1)=\left\{\begin{array}{cc}{r}_{i_1}\left(t\right)x_{i_1,j}\left(t\right)+r_{i_2}\left(t\right)x_{i_2,j}\left(t\right)-r_{i_3}\left(t\right)x_{i_3,j}\left(t\right)& j\&le;m\\ Great(AS,i,j)& j>m\end{array}\right.---\left(8\right)$

In formula: v_i,j(t+1) it is t+1 in period current pest population P_iFeature j state value；WithPoint Not Wei period t pest populationWithThe state value of feature j；i₁, i₂, i₃It is from { s₁, s₂..., s_LRandomly select in } Out, and meet i₁≠i₂≠i₃；Great (W, j) is meant that:

k₁, k₂, k₃It is the numbering of three different populations randomly selected from set W, i.e. meets k₁≠k₂≠k₃；

(S14) if 1/4 < q₀≤ 1/2, then perform natural enemy operator by formula (9), obtain v_i,j(t+1)；

$v_{i,j}(t+1)=\left\{\begin{array}{cc}{r}_i\left(t\right)x_{i,j}\left(t\right)+\underset{s\&Element;AS}{\max }\left\{r_s\left(t\right)x_{s,j}\left(t\right)\right\}-\underset{s\&Element;AS}{\max }\left\{r_s\left(t\right)x_{s,j}\left(t\right)\right\}& j\&le;m\\ Great(AS,i,j)& j>m\end{array}\right.---\left(9\right)$

(S15) if 1/2 < q₀≤ 3/4, then perform advantage operator as j≤m by formula (10), obtain v_i,j(t+1)；As j > m time by formula (11) perform advantage operator, obtain v_i,j(t+1)；

$v_{i,j}(t+1)=\left\{\begin{array}{cc}{r}_{g_a}\left(t\right)x_{g_a,j}\left(t\right)+r_{g_b}\left(t\right)x_{g_b,j}\left(t\right)-r_{g_c}\left(t\right)x_{g_c,j}\left(t\right)& \left|PM\right|\&GreaterEqual;3\\ r_i\left(t\right)x_{i,j}\left(t\right)+r_{g_a}\left(t\right)x_{g_a,j}\left(t\right)-r_{g_b}\left(t\right)x_{g_b,j}\left(t\right)& \left|PM\right|=2\\ r_{g_a}\left(t\right)x_{g_a,j}\left(t\right)& \left|PM\right|=1\\ x_{i,j}\left(t\right)& \left|PM\right|=0\end{array}\right.---\left(10\right)$

$v_{i,j}(t+1)=\left\{\begin{array}{cc}Great(PM,i,j)& \left|PM\right|\&GreaterEqual;1\\ x_{i,j}\left(t\right)& \left|PM\right|=0\end{array}\right.---\left(11\right)$

In formula, g_a、g_b、g_cAt { g₁, g₂..., g_LRandomly choose in }, it is desirable to g_a≠g_b≠g_c；

(S16) if 3/4 < q₀≤ 1, then as j≤m by the surging operator of formula (12), obtain v_i,j(t+1)；As j > m time hold by formula (13) The surging operator of row, obtains v_i,j(t+1)；

$v_{i,j}(t+1)=\left\{\begin{array}{cc}{r}_{h_a}\left(t\right)x_{h_a,j}\left(t\right)+r_{h_b}\left(t\right)x_{h_b,j}\left(t\right)-r_{h_c}\left(t\right)x_{h_c,j}\left(t\right)& \left|SM\right|\&GreaterEqual;3\\ r_i\left(t\right)x_{i,j}\left(t\right)+r_{h_a}\left(t\right)x_{h_a,j}\left(t\right)-r_{h_b}\left(t\right)x_{h_b,j}\left(t\right)& \left|SM\right|=2\\ r_{h_a}\left(t\right)x_{h_a,j}\left(t\right)& \left|SM\right|=1\\ x_{i,j}\left(t\right)& \left|SM\right|=0\end{array}\right.---\left(12\right)$

$v_{i,j}(t+1)=\left\{\begin{array}{cc}Great(SM,i,j)& \left|SM\right|\&GreaterEqual;1\\ x_{i,j}\left(t\right)& \left|SM\right|=0\end{array}\right.---\left(13\right)$

In formula, h_a、h_b、h_cAt { h₁, h₂..., h_LRandomly choose in }, it is desirable to h_a≠h_b≠h_c；

(S17) if p > E₀, then v is made_i,j(t+1)=x_i,j(t)；

(S18) make j=j+1, if j≤n, then go to step (S10), otherwise go to step (S19)；

(S19) perform accretive operatos by formula (14), obtain X_i(t+1)；

In formula:

X_i(t)=(x_i,1(t),x_i,2(t),…,x_i,n(t))；

V_i(t+1)=(v_i,1(t+1),v_i,2(t+1),…,v_i,n(t+1))；

(S20) make i=i+1, if i≤N, then go to step (S7), otherwise go to step (S21)；

(S21) if newly obtained globally optimal solution X^*t+1And the error between the globally optimal solution that the last time obtains meets minimum Require ε, then go to step (S24), otherwise go to step (S22)；

(S22) newly obtained globally optimal solution X is preserved^*t+1；

(S23) make t=t+1, if t≤G, then turn above-mentioned steps (S4), otherwise go to step (S24)；

(S24) terminate.