Abstract Discrete Optimization An alternative framework to Lagrangian relaxation

合集下载

Particle swarm approach for structural design optimization

Particle swarm approach for structural design optimizationR.E.Pereza,*,K.Behdinanba Institute for Aerospace Studies,University of Toronto,4925Duﬀerin Street,Toronto,Ont.,Canada M3H 5T6bDepartment of Aerospace Engineering,Ryerson University,350Victoria Street,Toronto,Ont.,Canada M5B 2K3Received 17February 2006;accepted 31October 2006Available online 6March 2007AbstractThis paper presents in detail the background and implementation of a particle swarm optimization algorithm suitable for constraint structural optimization tasks.Improvements,eﬀect of the diﬀerent setting parameters,and functionality of the algorithm are shown in the scope of classical structural optimization problems.The eﬀectiveness of the approach is illustrated by three benchmark structural optimization tasks.Results show the ability of the proposed methodology to ﬁnd better optimal solutions for structural optimization tasks than other optimization algorithms.Ó2007Elsevier Ltd.All rights reserved.Keywords:Structural optimization;Truss structures;Particle swarms;Constrained optimization;Soft Computing;Stochastic optimization1.IntroductionOver the past decade a number of optimization algo-rithms have been used extensively in structural optimiza-tion task,from gradient-based algorithms using convex and continuous design spaces,to non-gradient probabilis-tic-based search algorithms widely applied for global and non-convex design exploration.In the latter category many of these algorithms have been developed by mimicking nat-ural phenomena like simulated annealing [1],genetic algo-rithms [2],and bacterial foraging [3]among others.Recently,a family of optimization algorithms has been developed based on the simulation of social interactions among members of a speciﬁc species looking for food sources.From this family,the two most promising algo-rithms are ant colony optimization,which is based on the pheromone pathways used by ants to guide other ants in a colonies towards food sources [4],and particle swarm optimization or PSO,which is based on the social behav-iour reﬂected in ﬂock of birds,bees,and ﬁsh that adjusttheir physical movements to avoid predators,and to seek the best food sources [5].The PSO algorithm was ﬁrst proposed by Kennedy and Eberhart [6].It is based on the premise that social sharing of information among members of a species oﬀers and evo-lutionary advantage.A number of advantages with respect to other algorithms make PSO an ideal candidate to be used in optimization tasks.The algorithm is robust,well suited to handle non-linear,non-convex design spaces with discontinuities.It can handle continuous,discrete and inte-ger variable types with ease.As compared to other robust design optimization methods PSO is more eﬃcient,requir-ing fewer number of function evaluations,while leading to better or the same quality of results [7,8].In addition,it easiness of implementation makes it more attractive as it does not required speciﬁc domain knowledge information,internal transformation of variables or other manipula-tions to handle constraints.Furthermore,it is a popula-tion-based algorithm,so it can be eﬃciently parallelized to reduce the total computational eﬀort.Recently,the PSO has been proved useful on diverse engineering design applications such as logic circuit design [9],control design [10–12],and power systems design [13–15]among others.Applications in structures had been done in the area of structural shape optimization [16,17],0045-7949/$-see front matter Ó2007Elsevier Ltd.All rights reserved.doi:10.1016/pstruc.2006.10.013*Corresponding author.E-mail addresses:rperez@utias.utoronto.ca (R.E.Perez),kbehdina@acs.ryerson.ca (K.Behdinan)./locate/compstrucComputers and Structures 85(2007)1579–1588and in topology optimization[18],with promising results in such structural design applications.This paper focuses on the implementation and application of PSO for structural optimization.The general PSO methodology as well as the enhancement of the basic algorithm is introduced. Application of the algorithm to classical constraint struc-tural optimization problems is shown.The development of the paper is as follows:Section2presents the general formulation of the particle swarm optimization approach. In Section3analytical properties of the algorithm are discussed and improvement of the basic algorithms are pre-sented.Section4presents diﬀerent structural optimization case studies.They are used to analyze the behaviour and sensitivity of the PSO parameters,and demonstrate the eﬀectiveness of the approach inﬁnding optimal structural optimization solutions.Finally,the paper closes with some concluding remarks.2.Particle swarm algorithm2.1.Mathematical formulationThe particle swarm process is stochastic in nature;it makes use of a velocity vector to update the current posi-tion of each particle in the swarm.The velocity vector is updated based on the‘‘memory’’gained by each particle, conceptually resembling an autobiographical memory,as well as the knowledge gained by the swarm as a whole. Thus,the position of each particle in the swarm is updated based on the social behaviour of the swarm which adapts to its environment by returning to promising regions of the space previously discovered and searching for better positions over time.Numerically,the position x of a particle i at iteration k+1is updated as shown in Eq.(1)and illustrated in Fig.1.x i kþ1¼x ikþv ikþ1D tð1Þwhere v ikþ1is the corresponding updated velocity vector,and D t is the time step value[19].Throughout the present work a unit time step is used.The velocity vector of each particle is calculated as shown in Eq.(2),v i kþ1¼wv ikþc1r1ðp ikÀx ikÞD tþc2r2ðp gkÀx ikÞD tð2Þwhere v ikis the velocity vector at iteration k,r1and r2rep-resents random numbers between0and1;p ikrepresents thebest ever particle position of particle i,and p gkcorrespondsto the global best position in the swarm up to iteration k.The remaining terms are problem dependent parameters;for example,c1and c2represent‘‘trust’’parameters indicat-ing how much conﬁdence the current particle has in itself(c1or cognitive parameter)and how much conﬁdence ithas in the swarm(c2or social parameter),and w is the iner-tia weight.This later term plays an important role in thePSO convergence behaviour since it is employed to controlthe exploration abilities of the swarm.It directly impactsthe current velocity,which in turn is based on the previoushistory of rge inertia weights allow for widevelocity updates allowing to globally explore the designspace,while small inertia values concentrate the velocityupdates to nearby regions of the design space.2.2.Algorithm descriptionBased on the particle and velocity updates as explainedabove,the algorithm is constructed as follows:(1)Initialize a set of particles positions x iand velocitiesv irandomly distributed throughout the design spacebounded by speciﬁed limits.(2)Evaluate the objective function values fðx ikÞusing thedesign space positions x ik.(3)Update the optimum particle position p ikat currentiteration(k)and global optimum particle position p gk.(4)Update the position of each particle using its previousposition and updated velocity vector as speciﬁed inEqs.(1)and(2).(5)Repeat steps2–4until the stopping criteria ismet.For the current implementation the stopping cri-teria is deﬁned based on the number of iterationsreached.Eqs.(3)and(4)are used to obtain the initial position andvelocity vectors randomly distributed throughout thedesign space for each particle of the swarm.The term rin both equations represents a random number between0and1,while x min and x max represent the design variableslower and upper bounds respectively.x i¼x minþrðx maxÀx minÞð3Þv i¼x minþrðx maxÀx minÞD tð4Þ2.3.Algorithm analysis2.3.1.PSO search direction and stochastic step sizeLet us replace the PSO velocity update equation(2)intoEq.(1)to get the following expression:x ikþ1¼x ikþwV ikþc1r1ðp ikÀx ikÞD tþc2r2ðp gkÀx ikÞD tD tð5Þ1580R.E.Perez,K.Behdinan/Computers and Structures85(2007)1579–1588Factorizing the cognitive and social terms from the above equation we obtain the following general equation:x i k þ1¼x i k þwV i k D t þðc 1r 1þc 2r 2Þc 1r 1p i kþc 2r 2p g k c 1r 1þc 2r 2Àx i k ð6Þwhich follows the general gradient line-search form x i k þ1¼^x i k þa p k such that:^x i k ¼x i k þwV i k D t a ¼c 1r 1þc 2r 2 p k ¼c 1r 1p i kþc 2r 2p g k1122Àx i kð7Þwhere a can be regarded as a stochastic step size,and p k as a stochastic search direction.The stochastic step size will be limited only by the selection of social and cognitive param-eters;knowing that r 1,r 22[0,1]it will belong to the inter-val [0,c 1+c 2]with a mean value of ½c 1þc 2.Similarly,the stochastic search direction will belong to the interval½Àx i k ;c 1p i kþc 2p g k12Àx i k ,which is dependent not only on the so-cial and cognitive parameters,but also in the individual and global best positions.2.3.2.PSO convergenceBy replacing the particle velocity equation (2)into the position equation (1)and re-arranging the position term,the general form for the i th particle position equation at iteration k +1can be obtained asx i k þ1¼x i k ð1Àc 1r 1Àc 2r 2ÞþwV i k D t þc 1r 1p i k þc 2r 2p g kð8ÞSimilarly,re-arranging the position term in the particle velocity equation (2)leads to Vi k þ1¼Àx i kðc 1r 1þc 2r 2ÞD t þwV ik þc 1r 1p i k D t þc 2r 2p g kD tð9ÞEqs.(8)and (9)can be combined and written in matrixform asx i k þ1V i k þ1"#¼1Àc 1r 1Àc 2r 2w D t Àðc 1r 1þc 2r 2ÞD tw 2435x i k V i k þc 1r 1c 2r 2c 1r 1D t c 2r 2D t "#p i kp g k ð10Þwhich can be considered as a discrete-dynamic system rep-resentation for the PSO algorithm where [x i ,V i ]T is the state subject to an external input [p i ,p g ]T ,and the ﬁrst and second matrices correspond to the dynamic and input matrices,respectively.Assuming there is no external excitation in the dynamic system so [p i ,p g ]T is constant (other particles do not ﬁnd better positions),then a convergent behaviour could be maintained.If such is the case,as the iterations k !1then ½x i k þ1;V i k þ1 T ¼½x i k ;V i k T which reduces the dynamic metric form to0 ¼Àðc 1r 1þc 2r 2Þw D t Àðc 1r 1þc 2r 2ÞD tw À12435x i k V i k þc 1r 1c 2r 2c 1r 1D t c 2r 2D t "#p i k p g kð11Þwhich is true only when V i k ¼0and both x i k and p i k coincide with p g k .Note that such position is not necessarily a local or global minimizer,it is instead an equilibrium point.Such point,however,will improve towards the optimum if there is external excitation in the dynamic system,so better indi-vidual-best and global-best positions are found from the optimization process.2.3.3.PSO stabilityThe dynamic matrix eigenvalues derived from Eq.(10)can provide additional insight in the system stability and dynamic behaviour.The dynamic matrix characteristic equation is derived ask 2Àðw Àc 1r 1Àc 2r 2þ1Þk þw ¼0ð12Þwhere the eigenvalues are given ask 1;2¼1þw Àc 1r 1Àc 2r 2Æﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃð1þw Àc 1r 1Àc 2r 2Þ2À4w q 2ð13ÞThe necessary and suﬃcient condition for stability of a dis-crete-dynamic system is that all eigenvalues (k )derived from the dynamic matrix lies inside a unit circle around the origin on the complex plane,so j k i =1,...,n j <1.An analysis on the derived eigenvalues (12)leads to the follow-ing set of stability conditions for the PSO:c 1r 1þc 2r 2>0ðc 1r 1þc 2r 2Þ2Àw <1w <1ð14ÞFrom the above conditions the inertia term lies in the inter-val ðc 1r 1þc 2r 2ÞÀ1<w <1from which the upper limit for c 1r 1+c 2r 2can be derived as (c 1r 1+c 2r 2)<4.Knowing that r 1,r 22[0,1]the following parameter selection heuristic can be derived:0<ðc 1þc 2Þ<4ðc 1þc 2Þ2À1<w <1ð15ÞTherefore,if w ,c 1,and c 2are selected with the heuristics speciﬁed in (15),and meet the conditions speciﬁed in (14)the system has guaranteed convergence to an equilibrium point as discussed in Section 2.3.2.3.Algorithm improvements 3.1.Inertia weight updateDue to the importance of the inertia weight in control-ling the global/local search behaviour of the PSOR.E.Perez,K.Behdinan /Computers and Structures 85(2007)1579–15881581algorithm,a dynamic improvement has proven useful by forcing an initial global search with a high inertia weight (w%1)and subsequently narrowing down the algorithm exploration to feasible areas of the design space by decreas-ing its value towards local search values(w<0.5).Two main approaches had been used to deal with such dynamic update.In Shi and Eberhart[19],a dynamic variation of inertia weight is proposed by linearly decreasing w with each algorithmic iteration as shown in Eq.(16)w kþ1¼w maxÀw maxÀw mink maxkð16Þwhile in Fourie and Groenwold[16]a dynamic decrease of w value has been suggested based on a fraction multiplier (k w)as shown in Eq.(17),when no improvements had been made for a predeﬁned number of consecutive design iterationsw kþ1¼k w w kð17Þ3.2.Violated design points redirectionAn additional improvement introduced by Venter and Sobieszczanski-Sobieski[20]when dealing with violated constraint particles is included in the present optimization eﬀort.The method restricts the velocity vector of a violated particle to a usable feasible direction that would reduce the objective function while pointing back towards feasible regions of the design space.The modiﬁcation is made to the velocity vector of a particle i speciﬁed in Eq.(2)when there is one or more violated constraints by re-setting the velocity vector of the particle according to Eq.(18)which include only the obtained particle self-information of the best point and the current best point in the swarm.This strategy will introduce the new velocity vector pointing towards the feasible regions of the design space as shown in Fig.2.A new position for the violated constraint parti-cles is deﬁned by using Eq.(1)with the velocity vector modiﬁed asv i kþ1¼c1r1ðp ikÀx ikÞD tþc2r2ðp gkÀx ikÞD tð18ÞThe velocity of particle i at iteration k+1is thus only inﬂuenced by the best point found so far for that particle,and the current best point in the swarm.If both of these best points are feasible,the new velocity vector will point back to a feasible region of the design space.Otherwise, the new velocity vector will point to a region of the design space that resulted in smaller constraint violations.The re-sult is to have the violated particle move back towards the feasible region of the design space,or at least closer to its boundary,in the next design iteration.3.3.Constraints handlingSimilar to other stochastic optimization methods,the PSO algorithm is deﬁned for unconstrained problems.To accommodate the inclusion of constraints,a parameter-less adaptive penalty scheme is used as seen in Eq.(19).It uses the swarm information,such as the average of the objective function and the level of violation of each constraint during each iteration,in order to deﬁne diﬀerent penalties for dif-ferent constraints:f0ðxÞ¼fðxÞif x is feasible;fðxÞþP mi¼1k i g iðxÞotherwise8<:ð19ÞThe penalty parameter is deﬁned at each iteration ask i¼j fðxÞjg iðxÞPl¼1½ g lðxÞð20Þwhere f(x)is the objective function,m is the number of con-straints,g i(x)is a speciﬁc constraint value(with violated constraints having values larger than zero), fðxÞis the aver-age of the objective function values in the current swarm, and g iðxÞis the violation of the l th constraint averaged over the current population.The above formulation distributes the penalty coeﬃ-cients in a way that those constraints which are more diﬃ-cult to be satisﬁed will have a relatively higher penalty coeﬃcient.Such distribution is achieved by making the i th coeﬃcient proportional to the average of the violation of the i th constraint by the elements of the current popula-tion.An individual in the swarm whose i th violation equals the average of the i th violation in the current population for all i,has a penalty equal to the absolute value of the1582R.E.Perez,K.Behdinan/Computers and Structures85(2007)1579–1588average objective function of the population.Similarly,the average of the objective function equals fðxÞþj fðxÞj.4.Numerical examplesThe eﬀectiveness of the implemented PSO algorithm on structural optimization is shown through the use of four classical truss optimization examples.4.1.Six-node trussTheﬁrst example considers a well-known test problem corresponding to a10-bar truss non-convex optimization as describe by Sunar and Belegundu[21]and shown in Fig.3.In this problem the cross-sectional area for each of the10members in the structure are being optimized towards the minimization of total weight.The cross-sec-tional area varies between0.1in.2and35.0in.2.Con-straints are speciﬁed in terms of stress and displacement of the truss members.The allowable stress for each mem-ber is25,000psi for both tension and compression,and the allowable displacement on the nodes is±2in.,in the x and y directions.The density of the material is0.1lb/in.3, Youngs modulus is E=104ksi and vertical downward loads of100kips are applied at nodes2and4.In total the problem has a variable dimensionality of10and con-straint dimensionality of34(10tension constraints,10 compression constraints,and12displacement constraints). The eﬀect of the diﬀerent PSO parameters,speciﬁcally the eﬀect of social and cognitive parameters and diﬀerent dynamic inertia weight updates are presented,followed by general result of the sizing problem.4.1.1.Social and cognitive parameters variation eﬀectA more detailed analysis of the eﬀect of c1and c2param-eters is necessary to determine the sensitivity of the param-eters in the overall optimization procedure.The eﬀect of the cognitive(c1)and social(c2)parameters on the general optimization history can be seen in Fig.4where objective function values vs.iteration number is shown.The results are representative of more than20trials for each tested case,with aﬁxed inertia weight value of0.875,and the same initial particles and velocity values.The parameters were varied according to the formula c2=4Àc1,where c1was varied in the interval[4,0].It can be seen clearly that when only cognitive(c1=4,c2=0)or only social values (c1=0,c2=4)are used the resulting history converges very fast,within theﬁrst10iterations,but do not improve after the initial convergence indicating that the algorithm reaches a local suboptimal and is not able to escape from it,due to the lack of information exchange either from social or cognitive sources respectively.When higher emphasis is placed on the social exchange of information (c1=1,c2=3)the algorithm is able to gradually converge to better regions of the design space,but again if a local optimum is found it cannot escape from it.This fact indi-cates an overshooting of the algorithm where particle posi-tions are updated based mostly on social information not promoting the exchange of individual information which could be useful when an individual particle has found a bet-ter region of the design space.When the cognitive and social parameters are in balance,(c1=c2=2,so the mean of Eq.(2)is equal to1)or with slightly higher cognitive value(c1=3,c2=1),solutions converge near the global optimum solution;these two situations are shown by the gradual decrease of objective function values and subse-quent improvements over time.Notice that a slightly higher value for the cognitive parameter results in better solutions.This result is due to the fact that individuals con-centrate more in their own search regions avoiding over-shooting the best design space region,but promoting some global exchange of information to make the swarm to point toward the best regions.4.1.2.Inertia weight variation eﬀectTo analyze the eﬀect of the inertia weight variation, three diﬀerent test cases were analyzed.Theﬁrst one con-siders aﬁxed mean inertia weight of the typical w interval [1.00,0.35],the second one uses a linearly varying update over the same interval,while the third one uses a dynamic update as described by Eq.(4)with a fraction multiplier k w=0.975.Results are given in Fig.5.Notice thatﬁxedR.E.Perez,K.Behdinan/Computers and Structures85(2007)1579–15881583inertia weight allow for a fast convergence(within50iter-ations)to sub-optimal regions of the design space while lin-ear and dynamic weight updates converge slowly,due to the initial global search nature force by higher weight update values,and gradually converges to optimal regions of the design space approximately after iteration120;no clear distinction in the solution is apparent between both varying inertia weight updates,but the dynamic update tend to converge faster to the optimal solution which could be useful if a speciﬁc converge value is prescribed as the stop criteria of the algorithm.It is important to note that the optimal selection of the PSO parameters are in general problem-dependent.How-ever,the obtained results conﬁrms the expectations as derived from the theoretical analysis(see i.e.[22–25])and experiments(see i.e.[26–28])regarding the sensitivity and behaviour of such tuning parameters.As long the stability conditions presented in(14)are met,it is observed that maintain an approximately equal weighting between the cognitive and social parameters in the interval of1.4–2lead to the optimal convergence for the PSO.Similarly,aﬁxed inertia weight makes the PSO to stagnated with a high rate of convergence in local sub-optimal values,while dynami-cally decreasing such parameter allows for better explora-tion of the design space with similar rates of convergence.4.1.3.Structural optimization resultsTable1shows the PSO best and worst results of20inde-pendent runs for the10-bar truss problem.Based on the above inertia study,results were obtained using a dynamic inertia weight update to provide the best global/local opti-mization trade.Other published results,for the same case, are summarized as well in Tables1and2.They are foundTable1Optimization results for the10-bar trussTruss area PSO best PSO worst Gellatly andBerke[30]Schimit andMiura[29]Ghasemi et al.[37]Schimit andFarshi[34]Dobbs andNelson[31]0133.50033.50031.35030.57025.73033.43230.500 020.1000.1000.1000.3690.1090.1000.100 0322.76633.50020.03023.97024.85024.26023.290 0414.41713.30415.60014.73016.35014.26015.428 050.1000.1000.1400.1000.1060.1000.100 060.1000.1000.2400.3640.1090.1000.210 077.534 6.82638.3508.5478.7008.3887.649 0820.46718.93522.21021.11021.41020.74020.980 0920.39218.81422.06020.77022.30019.69021.818 100.1000.1000.1000.3200.1220.1000.100 Weight5024.215176.275112.005107.305095.655089.005080.00 Table2Optimization results for the10-bar truss(continuation)Truss area Rizzi[32]Haug andArora[39]Haftka andGrdal[40]Adeli andKamal[35]El-Sayed andJang[36]Galante[38]Memari andFuladgar[33]0130.73130.03130.52031.2832.96630.44030.561 020.1000.1000.1000.100.1000.1000.100 0323.93423.27423.20024.6522.79921.79027.946 0414.73315.28615.22015.3914.14614.26013.619 050.1000.1000.1000.100.1000.1000.100 060.1000.5570.5510.100.7390.4510.100 078.5427.4687.4577.90 6.3817.6287.907 0820.95421.19821.04021.5320.91221.63019.345 0921.83621.61821.53019.0720.97821.36019.273 100.1000.1000.1000.100.1000.1000.100 Weight5061.605060.9205060.805052.005013.244987.004981.1 1584R.E.Perez,K.Behdinan/Computers and Structures85(2007)1579–1588using diﬀerent optimization approaches including gradient based algorithms both unconstrained [29],and constrain-ted [30–33],structural approximation algorithms [34],con-vex programming [35],non-linear goal programming [36],and genetic algorithms (global search optimization)[37,38].We can see that the PSO provide good results as compared with other methods for this problem.The opti-mal solutions found by the PSO meet all constraints and have only two active constraints including the displace-ments at nodes 3and 6.4.2.25-bar trussThe second example considers the weight minimization of a 25-bar transmission tower as described by Schmit and Fleury [41]and shown in Fig.6.The design variables are the cross-sectional area for the truss members,which are linked in eight member groups as shown in Table 3.Loading of the structure is presented in Table 4.Con-straints are imposed on the minimum cross-sectional area of each truss (0.01in.2),allowable displacement at each node (±0.35in.),and allowable stresses for the members in the interval [À40,40]ksi.In total this problem has a var-iable dimensionality of 8and a constraint dimensionality of 84.Table 5present the best and worst result of 20indepen-dent runs of the PSO after 200iterations with a swarm of 40individuals for the 25-bar truss optimization task.This table presents as well results for the same optimization task from diﬀerent research eﬀorts.Clearly,the PSO yield excel-lent solutions for both its best and worst results.The best PSO solution provides the lowest structural weight of any reported solution,while the worst PSO solution is in agree-ment with those reported in the literature.The optimal solutions obtained by the PSO have the same active con-straints as reported in other references as follows:the dis-placements at nodes 3and 6in the Y -direction for both load cases and the compressive stresses in members 19and 20for the second load case.Also,there are no con-straint violations.4.3.72-bar trussThe third example comprises the optimization of a four-story 72-bar space truss shown in Fig.7.The dimen-sions are indicated as well in the ﬁgure.The optimizationTable 3Truss members area grouping for the 25-bar truss Group Truss members A11A22–5A36–9A410,11A512,13A614–17A718–21A822–25Table 4Nodal loads for the 25-bar truss Node F x F y F z 11000À10,000À10,00020À10,000À10,0003500006600Table 5Optimization results for the 25-bar truss Area group PSO best PSO worst Zhou and Rozvany [42]Haftka and Grdal [40]Erbatur et al.[43]Zhu [44]Wu and Chow [45]A10.10000.10000.0100.0100.10.10.1A2 1.02270.9895 1.987 1.987 1.2 1.90.5A3 3.4000 3.4000 2.994 2.991 3.2 2.6 3.4A40.10000.10000.0100.0100.10.10.1A50.1000 3.40000.0100.012 1.10.1 1.5A60.63990.69990.6840.6830.90.80.9A7 2.0424 1.9136 1.677 1.6790.4 2.10.6A8 3.4000 3.4000 2.662 2.664 3.4 2.6 3.4Weight485.33534.84545.16545.22493.80562.93486.29R.E.Perez,K.Behdinan /Computers and Structures 85(2007)1579–15881585。

Abstract Discrete Optimization A hybrid genetic algorithm for the Three-Index Assignment Pr

Discrete OptimizationA hybrid genetic algorithm for the Three-IndexAssignment ProblemGaofeng Huang,Andrew Lim*Department of Industrial Engineering and Engineering Management,Hong Kong University of Science and Technology,Clear Water Bay,Kowloon,Hong KongReceived 18December 2003;accepted 29September 2004AbstractThe Three-Index Assignment Problem (AP3)is well-known problem which has been shown to be NP -hard.This problem has been studied extensively,and many exact and heuristic methods have been proposed to solve it.Inspired by the classical assignment problem,we propose a new local search heuristic which solves the problem by simplifying it to the classical assignment problem.We further hybridize our heuristic with the genetic algorithm (GA).Extensive experimental results indicate that our hybrid method is superior to all previous heuristic methods including those pro-posed by Balas and Saltzman [Operations Research 39(1991)150–161],Crama and Spieksma [European Journal of Operational Research 60(1992)273–279],Burkard et al.[Discrete Applied Mathematics 65(1996)123–169],and Aiex et al.[GRASP with path relinking for the three-index assignment problem,Technical report,INFORMS Journal on Computing,in press.Available from:</~mgcr/doc/g3index.pdf>].Ó2004Elsevier B.V.All rights reserved.Keywords:Combinatorial optimization;Heuristic;Genetic algorithms;Assignment1.IntroductionThe Three-Index Assignment Problem (AP3),also known as the 3-Dimensional Assignment Problem,has been a subject of extensive research since the sixties [10,11].The 0–1programming model for AP3isminXi 2I ;j 2J ;k 2Kc ijk x ijksubject to Xj 2J ;k 2Kx ijk ¼1;8i 2I ;ð1ÞXi 2I ;j 2Jx ijk ¼1;8k 2K ;ð2ÞXi 2I ;k 2Kx ijk ¼1;8j 2J ;ð3Þx ijk 2f 0;1g ;8f i ;j ;k g 2I ÂJ ÂK ;ð4Þ0377-2217/$-see front matter Ó2004Elsevier B.V.All rights reserved.doi:10.1016/j.ejor.2004.09.042*Corresponding author.E-mail addresses:gfhuang@ust.hk (G.Huang),iealim@ust.hk (A.Lim).European Journal of Operational Research xxx (2004)xxx–xxxwhere I ,J ,K are three disjoint sets with j I j =j J j =j K j =n .In fact,AP3can be considered as an optimiza-tion problem on a complete tripartite graph K n ,n ,n =(I [J [K ;(I [J )·(I [K )·(J [K )),where I ,J and K are the three disjoint vertex sets.The cost of choosing triangle (i ,j ,k )2I ·J ·K is c i ,j ,k .The objective of AP3is to choose n disjoint trian-gles (i ,j ,k )so that the total cost is minimized.As a result,of all the n 3entries of c i ,j ,k ,only n entries are chosen.If we line up these n chosen triples,we get three permutations of 1,2,...,n .This means that a solution to AP3can also be represented by three permutations.The permutation based AP3formu-lation is minX n i ¼1c I ði Þ;p ði Þ;q ði Þwith I ;p ;q 2p N ;where p N denotes the set of all permutations on the set of integers N ={1,2,...,n }.Since we do not care about the order of these n triples,we can reorder them such that I (i )=i ,or I =(1,2,3,...,n ).Once this ‘‘index permutation’’is ﬁxed,the solution to AP3can be represented by using a pair of permutations (p ,q ).AP3can be formulated as:Instance :a matrix C ={c i ,j ,k }n ·n ·n Solution :(p ,q ),p ,q 2p NObjective :to minimize C ðp ;q Þ¼Pn i ¼1c i ;p ði Þ;q ði ÞIt is quite obvious that AP3is a straightforward extension of the classical two-dimensional assign-ment problem (AP2)deﬁned below:Instance :matrix D ={d i ,j }n ·n (bipartite graph)Solution :q =(q 1,q 2,...,q n ),q 2p N Objective :to minimize D ðq Þ¼P n i ¼1d i ;q ði ÞAlthough it is well-known that AP2can be solved eﬃciently in polynomial time [8,9],the AP3is NP -hard since the 3-D Matching Problem,which is one of the basic NP -hard problems [7]is one of the special case of AP3.Both exact and heuristic algorithms have been proposed to solve AP3[1–6,10,11].Among these,Balas and Saltzman [5]proposed the MAX-REGRET and VARIABLE DEPTH INTER-CHANGE heuristics.Crama and Spieksma [4]studied a special case of AP3by restricting the cost of edges in any triangle to obey the triangle inequalities.Burkard et al.[2]focused on AP3with decomposable cost coeﬃcients.However,even for these two special cases,AP3is still NP -hard.Re-cently,Aiex et al.[1]applied GRASP with Path Relinking for AP3and obtained better results than all other existing heuristics.The purpose of this paper is to present a hybrid genetic algorithm for AP3.Genetic algorithm (GA)is one of the most successful evolutionary algorithms and is based on an analogy with Dar-win Õs evolution theory of natural selection and ‘‘survival of the ﬁttest’’.Genetic algorithms have shown to be a competitive technique to solve gen-eral combinatorial optimization problems [12].If we can incorporate problem-speciﬁc knowledge into GA,the results can be further improved.The hybridization between local search (LS)and GA reﬂects this idea and becomes popular recently [13].In this paper,we ﬁrst propose a new local search heuristic algorithm for AP3,which utilizes problem-speciﬁc knowledge to solve AP3.We then hybridize this local search heuristic with GA.Experiments indicate that GA beneﬁted from this hybridization.We test our algorithm on three clas-ses of standard benchmarks and report the compu-tational results.The rest of this paper is organized as follows.Section 2ﬁrst describes the details of our proposed local search heuristic algorithm.In Section 3,we hybridize this local search with genetic algorithm.The results of extensive experiments are reported in Section 4.Finally,we present our conclusions in Section 5.2.Local search 2.1.Project 3D onto 2DAlthough AP3is NP -hard,it is well-known that AP2,the classical two-dimensional assign-ment problem,can be solved in polynomial time.A straightforward min-cost network ﬂow imple-2G.Huang,A.Lim /European Journal of Operational Research xxx (2004)xxx–xxxmentation for AP2takes O (n 4)time,an eﬃcient implementation of the Hungarian algorithm can solve AP2in O (n 3)time [8,9].As stated earlier,any solution to AP3consists of two permutations p and q ,while a solution to AP2consists of only one permutation q .Given an initial solution (p ,q )for AP3,let d i ;j ¼c i ;p ði Þ;j ;8i ;j 21;2;...;n ;ð5Þwe getminp ;q 2p NX n i ¼1c i ;p ði Þ;q ði Þ¼min q 2p NX n i ¼1d i ;q ði Þ:Consequently,if we ﬁx permutation p ,the optimi-zation of q becomes an AP2problem,and vice versa.Therefore,our idea is to optimize one permuta-tion subject to the other permutation being ﬁxed.To illustrate this,we use an example (instance bs _4_5.dat from Balas and Saltzman Dataset,see Section 4.2.1).As shown in Fig.1,a random initial assignment costs 177.Fig.2shows the optimiza-tion of q by applying the Hungarian Method.The objective function decreases from 177to 72(the dotted lines show the original assignment,while the new assignment is shown by bold lines).2.2.Iterative local searchAs shown in Fig.2,we construct a bipartitegraph based on Eq.(5).Symmetrically,there are two other ways to construct a bipartite graph from a tripartite graph:let d i ;j ¼c i ;j ;q ði Þ;8i ;j 21;2;...;nð6Þorlet d i ;j ¼c j ;p ði Þ;q ði Þ;8i ;j 21;2;...;n :ð7ÞFig.3illustrates Eq.(6),while Fig.4corresponds to Eq.(7).Using the above,the three parts of one solution to AP3,namely permutation p ,permutation q and the index permutation I ,can be optimized sepa-rately and eﬃciently.As shown in Fig.5,our local search heuristic iteratively optimizes the three parts (shaded parts)until no more improvement can be achieved.Algorithm 1illustrates the implementation.Algorithm 1.Local search (LS)for AP31:(p ,q ) any initial solution 2:ﬁnishFlag false 3:while not ﬁnishFlagdoG.Huang,A.Lim /European Journal of Operational Research xxx (2004)xxx–xxx34:ﬁnishFlag true5:for part of Solution1to3do6:construct corresponding bipartite graph according Eqs.(5)–(7)7:optimize it by applying the Hungarian Method on the bipartite graph8:if Obj.Value decreases9:ﬁnishFlag false10:end if11:end for12:end while3.Hybrid genetic algorithmHolland[12]ﬁrst introduced the concept of ge-netic algorithms in the year1975.Genetic algo-rithms have shown to be a competitive technique for solving general combinatorial optimization problems,and many variations have been pro-posed during last three decades.Ref.[14]provides a review of GA for reference.Fig.6illustrates the structure of GA we used in this paper.The initial‘‘generation’’is randomly generated.The‘‘Crossover Operator’’randomly chooses two individuals(as parents)according to theirﬁtness function values,and generates one new individual.This newly generated individual will‘‘mutate’’with a very small probability.All of these newly generated individuals are put into a candidates pool,whose size is normally twice the parent population.The algorithm then applies the‘‘survival of theﬁttest’’principle to the candi-dates pool and selects the top half of the individu-als with goodﬁtness value from the candidate pool will be selected to form the next‘‘generation’’. After several generations,the GA terminates after satisfying some termination criteria.Based on this structure,we present the details of the components of the implementation of our ge-netic algorithm below:•Individual representation andﬁtness function: In AP3,any pair of permutations p=(p1p2...p n)and q=(q1q2...q n)is a feasiblesolution.It is very common to use a2n lengthvector as the representation(chromosome)ofan individual.However,as stated earlier,once we get permutation p,we can computepermutation q by the Hungarian Method.Therefore,in our GA,we only use permuta-tion p(a n length vector)as the chromosomeencoding.The objective function C(p,q)actsas theﬁtness function of chromosome.•Initial generation:For the initial generation, each individual is set to be a pair of ran-domly generated permutations,after whichthe proposed local search algorithm isapplied to improve its quality.Therefore,every individual in the initial generation isalready a local minimum in the solutionspace.•Termination criteria:Two termination crite-ria are used in our implementation.Theﬁrstis convergence detection.When selectingnext generation from the candidates pool,we use a duplicate detection scheme to guar-antee that the individuals we selected are alldiﬀerent.If the number of diﬀerent individu-als in the candidates pool is less thanneeded4G.Huang,A.Lim/European Journal of Operational Research xxx(2004)xxx–xxxpopulation,it means that the evolution pro-cess has converged.When this happens,thealgorithm terminates.The second termina-tion criterion is improvement detection.Ifthe solutionÕs quality does not improve in10generations,the algorithm terminates.•Crossover operator:Many crossover opera-tors have been proposed[14],includingPMX,CX and OX.However,there is no evi-dence that one crossover operator is superiorto another operator.We used only PMX inour algorithm.The Partially Mapped Crossover(PMX)opera-tor was proposed by Goldberg and Lingle[15]. The basic idea of PMX is to exchange a partial segment between two parents.An adjustment step is needed to maintain the feasibility of the new children(see Fig.7).•Mutation operator vs.hybridization with local search:We did not implement mutationoperators in our algorithm because it wasnot very useful.Others have reported similarexperiences with mutation in GA[16,17].Instead,we replace mutation operator byour local search module.When a new indi-vidual is generated from Crossover Operator,we apply local search to improve its qualitybefore putting it into the candidates pool.By this method,we hybridize local searchwith genetic algorithm.Fig.8illustrates thestructure of our hybridization(LSGA).4.Experiments4.1.Preliminary experimentsPreliminary experiments were conducted to tune our algorithm because the parameter setting can inﬂuence the performance of a genetic algo-rithm substantially.As a compromise between solutionÕs quality and running time,we set the population size to be100and the candidate pool size to be200.Fig.9shows the comparison between pure GA, a multi-round local search(LS),and LSGA on in-stance bs_26_1.dat from the Balas and Saltzman Dataset(see Section4.2.1).The following remarks can be made from Fig.9:G.Huang,A.Lim/European Journal of Operational Research xxx(2004)xxx–xxx5•Pure GA is obviously bad in performance.•multiLS is a multi-round local search algo-rithm.In each round,it starts with a randominitial permutation and uses local search toimprove the solution.multiLS has the abilitytoﬁnd relatively good solutions in a shorttime.This reﬂects the eﬀectiveness of ourlocal search algorithm.Unfortunately,evenif much more time is given,this algorithmcannot improve the best solution.For theinstance bs_26_1.dat,no better solutionscan be found after1second.•LSGA shows the power of hybridization of local search and genetic algorithm.It is capa-ble ofﬁnding very competitive solutions.Hence,it is evident that our hybridization of GA and LS is successful.Guided by the GA,it is possible for local search to improve the quality of solution consistently;and with the help of LS, genetic algorithm becomes competitive inﬁnding good solutions.putational resultsIn this section,we demonstrate the eﬀective-ness of our hybrid genetic algorithm(LSGA) by testing our heuristic on three benchmark data-sets.All the codes are implemented by C/C++un-der a Pentium III800MHz PC with128M memory.For each instance,our LSGA is run only once.However,the computing machines used in AiexÕs paper[1]is a SGI Challenge R10000ma-chine.Therefore,in order to compare the CPU time,a scaling scheme is used according to SPEC (Standard Performance Evaluation Corporation /osg/cpu2000/),which points out that PIII800is not more4times faster than SGI Challenge R10000(see AppendixA).4.2.1.Balas and Saltzman DatasetThis dataset is generated by Balas and Saltzman [5].It includes60test instances with the problem size n=4,6,8,...,22,24,26.For each n,ﬁve in-stances are randomly generated with the integer cost coeﬃcients c i,j,k uniformly distributed in the interval[0,100].Table1shows the results of our experiments on this dataset.Each row stores the average score of theﬁve instances with the same size.The Column ‘‘Optimal’’shows the optimal solution reported by Balas and Saltzman,while column‘‘B–S’’is the result of their VARIABLE DEPTH INTER-CHANGE heuristic.Column‘‘GRASP with Path Relinking’’is the result reported in[1].ColumnTable1Balas and Saltzman Dataset(12*5instances)n Optima B–S GRASP with Path Relinking multiLS LSGAAvg.obj. value Avg.obj.valueAvg.obj.valueAverage CPU time(seconds)Avg.obj.valueAvg.time(seconds)Avg.obj.valueAvg.time(seconds)R10000PIII800PIII800PIII800442.243.2–––42.20.0042.20.00 640.245.4–––40.20.0140.20.01 823.833.6–––33.60.0123.80.03 101940.8–––22.60.01190.37 1215.62415.674.79>18.726.20.0215.60.87 141022.410106.55>26.6426.40.0210 1.73 16102510.2143.89>35.9726.00.0310 1.89 18 6.417.67.4190.88>47.7224.60.037.2 2.95 20 4.827.4 6.4246.70>61.6826.80.04 5.2 4.01 22418.87.8309.64>77.4126.40.05 5.6 4.54 24 1.8147.4382.45>95.6123.20.06 3.2 5.66 26 1.315.78.4465.20>116.323.20.07 3.610.78 6G.Huang,A.Lim/European Journal of Operational Research xxx(2004)xxx–xxx‘‘multiLS’’is the result of our multi-round local search algorithm,which is terminated after100 rounds.Finally,Column‘‘LSGA’’shows the re-sult of our hybrid genetic algorithm.The best re-sults among these algorithms are underlined in the table.Our LSGA outperformed other heuristics. However,these instances are hard since there still exist visible gaps between LSGAÕs solutions and the optimal solutions for n P18.It is very interesting that the overall objective function value decreases with the growth of prob-lem size n.One possible explanation is as there are n3entries of coeﬃcients c i,j,k uniformly distributed in[0,100].This means,on average,there are n3‘‘0’’s,n3‘‘1’’s,etc.On the other hand,only n en-tries of coeﬃcients c i,j,k will be chosen and counted in the objective function value.Therefore,with the growth n,it becomes easier toﬁnd the assignment with low cost.It is evident that our LSGA can provide much better solutions than GRASP and B–S.Further-more LSGA is about10times faster than GRASP in these instances.4.2.2.Crama and Spieksma DatasetCrama and Spieksma generated this dataset by restricting coeﬃcients c i,j,k=d i,j+d i,k+d j,k[4]. There are three types of instances in this dataset. For each type,three instances of size n=33and three instances of size n=66are generated.Tables2–4report the experimental results of the three types,with a total of18instances.In these tables,column‘‘C–S’’reports the result of Crama and SpieksmaÕs heuristic.For all these18instances,our LSGA can al-ways provide the best solutions over all heuristics. Obviously,these solutions are very close to opti-mal since the deviation from the Lower Bound is rather small.Table2Crama and Spieksma Dataset,Type I(six instances)CaselD n Lower Bound C–S GRASP with Path Relinking multiLS LSGAAvg.obj. value Avg.obj.valueAvg.obj.valueAverage CPU time(seconds)Avg.obj.valueAvg.time(seconds)Avg.obj.valueAvg.time(seconds)R10000PIII800PIII800PIII8003DA99N133160716181608660.5>165.1316080.0316080.03 3DA99N233139514111401680.5>170.1314010.0214010.11 3DA99N333160416091604676.1>169.0316040.0316040.11 3DA198N16626542668266415470.1>3867.526620.2126620.55 3DA198N26624332469244915010.9>3752.724490.2124490.27 3DA198N36627482775275915084.6>3771.227580.2227580.58 Table3Crama and Spieksma Dataset,Type II(six instances)CaselD n Lower Bound C–S GRASP with Path Relinking multiLS LSGAAvg.obj. value Avg.obj.valueAvg.obj.valueAverage CPU time(seconds)Avg.obj.valueAvg.time(seconds)Avg.obj.valueAvg.time(seconds)R10000PIII800PIII800PIII8003DIJS9N133477248614797766.06>191.5247970.0647970.113DIJ99N233503551425068772.84>193.2150680.1350670.263DIJ99N333426043524287762.19>190.5542870.0842870.263DI198N16696339780969414629.1>3657.396870.519684 4.863DI198N26688319142895414922.9>3730.789470.538944 3.353DI198N36696709888975114391.7>3597.997470.539745 3.09G.Huang,A.Lim/European Journal of Operational Research xxx(2004)xxx–xxx7It is surprising that for these instances,LSGA is about1000times faster than GRASP.GRASP needs several hours to get the solutions while LSGA only takes several seconds.4.2.3.Burkard,Rudolf and Woeginger DatasetBurkard et al.[2]described this dataset with decomposable cost coeﬃcients,which means that c i,j,k=a iÆb jÆc k.For each problem size n=4,6, 8,...,16,100test instances are provided.Table5is the result statistics of all these700in-stances,where each row is the average of the100 instances with same size n;column‘‘B–R–W’’re-ports the result of Burkard et al.Õs heuristic.For these test instances,our LSGA provides the same results with GRASP with about100times faster speed.However,this dataset is considered to be easy since even multiLS can also oﬀer very competitive solutions,which is even faster.5.ConclusionsThe Three-Index Assignment Problem(AP3)is studied in this paper.We proposed a new local search heuristic for AP3and hybridized it with ge-netic algorithm.Experiments indicated that this hybridization is successful.Computational results showed that our hybrid genetic algorithm(LSGA)outperformed other existing heuristics.For those benchmark instances that we have tested,LSGA is about10–1000times faster than GRASP and can always oﬀer the best solutions in several seconds.Appendix AAccording to Table6,PIII800:R10000=(386/ 233)*(20.7/8.75)=3.91<4.Table4Crama and Spieksma Dataset,Type III(six instances)CaselD n Lower Bound C–S GRASP with Path Relinking multiLS LSGAAvg.obj. value Avg.obj.valueAvg.obj.valueAverage CPU time(seconds)Avg.obj.valueAvg.time(seconds)Avg.obj.valueAvg.time(seconds)R10000PIII800PIII800PIII8003D1299N133133135133490.79>122.71330.011330.01 3D1299N233130137131471.21>117.81320.011310.03 3D1299N333130135131451.72>112.931310.011310.02 3D1198N1662832932865322.97>1330.72870.052860.15 3D1198N2662812942865126.86>1281.72860.052860.16 3D1198N3662802932825059.06>1264.82830.052820.23 Table5Burkard,Rudolf and Woeginger Dataset(700instances)n B–R–W GRASP with Path Relinking multiLS LSGAAvg.obj. value Avg.obj.valueAverage CPU time(seconds)Avg.obj.valueAvg.time(seconds)Avg.obj.valueAvg.time(seconds) R10000PIII800PIII600PIII8004443.7–––433.60.00443.60.00 6634.2–––633.720.00633.720.01 8819.94–––819.160.01819.160.03 10960.55–––959.420.03959.410.07 121188.021186.8168.3>17.11186.830.041186.810.13 141469.271467.7498.1>24.51467.760.071467.740.23 161476.991475.13139.3>34.81475.150.101475.130.40 8G.Huang,A.Lim/European Journal of Operational Research xxx(2004)xxx–xxxReferences[1]R.M.Aiex,M.G.C.Resende,P.M.Pardalos,G.Toraldo,GRASP with path relinking for the three-index assignment problem,Technical report,INFORMS Journal on Com-puting,in press.Available from:<http://www.research.att.com/~mgcr/doc/g3index.pdf>.[2]R.E.Burkard,R.Rudolf,G.J.Woeginger,Three dimen-sional axial assignment problems with decomposable cost coeﬃcients,Discrete Applied Mathematics65(1996)123–169.[3]R.E.Burkard,R.Rudolf,Computational investigationson3-dimensional axial assignment problems,Belgian Journal of Operations Research32(1993)85–98.[4]Y.Crama,F.C.R.Spieksma,Approximation algorithmsfor three-dimensional assignment problems with triangle inequalities,European Journal of Operational Research60 (1992)273–279.[5]E.Balas,M.J.Saltzman,An algorithm for the three-indexassignment problem,Operations Research39(1991)150–161.[6]P.Hansen,L.Kaufman,A primal–dual algorithm for thethree-dimensional assignment problem,Cahiers du CERO 15(1973)327–336.[7]M.R.Garey,D.S.Johnson,Computers and intractabil-ity—A guide to the theory of NP-completeness,W.H.Freeman and Company,1979.[8]H.W.Kuhn,The Hungarian Method for the assignmentproblem,Naval Research Logistics Quarterly2(1955)83–97.[9]wler,Combinatorial Optimization:Networks andMatroids,Holt,Reinhart and Winston,New York,1976.[10]W.P.Pierskalla,The tri-substitution method for the three-dimensional assignment problem,CORS Journal5(1967) 71–81.[11]W.P.Pierskalla,The multidimensional assignment prob-lem,Operations Research16(1968)422–431.[12]J.Holland,Adaptation in Natural and Artiﬁcial Systems,University of Michigan Press,Ann Arbor,MI,1975. [13]W.E.Hart,R.K.Belew,Optimization with genetic algo-rithms hybrids that use local search,in:R.K.Belew,M.Mitchell(Eds.),Adaptive Individuals in Evolving Popula-tions,Addison-Wesley,Reading,MA,1996.[14]rran˜aga,C.M.H.Kuijpers,R.H.Murga,I.Inza,S.Dizdarevic,Genetic algorithms for the travelling salesman problem:A review of representations and operators, Artiﬁcial Intelligence Review13(1999)129–170.[15]D.Goldberg,R.Lingle,Alleles,loci,and the travelingsalesman problem,in:Proceedings of the1st International Conference on Genetic Algorithms and their Applications, Lawrence Erlbaum,Hillsdale,NJ,1985,pp.154–159. [16]D.Goldberg,Genetic Algorithms in Search,Optimizationand Machine Learning,Addison-Wesley,Reading,MA, 1989.[17]G.Syswerda,Schedule optimization using genetic algo-rithms,in:L.Davis(Ed.),Handbook of Genetic Algo-rithms,Van Nostrand Reinhold,New York,1991,pp.332–349.Table6CPU Benchmark from SPEC(/osg/ cpu2000/)Intel PIII800IntelPIII500SGI ChallengeR10000Benchmark CPU9520.78.75Benchmark CPU2000386233G.Huang,A.Lim/European Journal of Operational Research xxx(2004)xxx–xxx9。

【doc】有限差分法在一维输运方程定解中的运用

有限差分法在一维输运方程定解中的运用2012年2月第12卷第1期廊坊师范学院(自然科学版)JournalofLangfangTeachersCollege(NaturalScienceEdition)Feb.2012V0l_12No.1有限差分法在一维输运方程定解中的运用林喜季(福建江夏学院,福建福州350108)【摘要】有限差分方法就是一种数值解法,在一维输运方程定解中可以巧用它来解题,把表示变量连续变化关系的偏微分方程离散为有限个代数方程,然后利用电子计算机求此线性代数方程组的解.【关键词】一维输运方程;有限差分法;定解问题FiniteDifferenceMethodin0ne—DimensionalTransport EquationintheUseofDefiniteSolutionLINXii【Abstract】Thefinitedifferencemethodisanumericalmethod,one-dimensionaltransportequationinth esolutioncanbeskillfullyusedittosolveproblems,thecontinuousvariationofthatvariablepartialdifferential equationisdiscretizedintoafinitenumberofalgebraicequations,andthenusethecomputertosolveThishnearalgebraice quations.【Keywords】one-dimensionaltransportequation;finitedifferencemethod;definitesolutionproblem[中图分类号]0175.2(文献标识码]A[文章编号]1674—3229(2012)01—0016—03 在数学中,有限差分法的内涵是指用差商代替微商,即用泰勒级数展开式将变量的导数写成变量在不同时间或空间点值的差分形式的方法.它的基本思想是按时间步长和空间步长将时间和空间区域剖分成若干方格网,用泰勒级数展开近似式代替所用偏微分方程中出现的各阶导数,从而把表示变量连续变化关系的偏微分方程离散为有限个代数方程,然后,解此线性代数方程组.l导数用泰勒级数展开近似式导数(微商)y:=m0=±,是无限小的微分m0△),除以无限小的微分是△的商.它可以分别近似为:,,=dxAx=(1)y=Ax=(2)y=dxAx=坐(3)式(1),(2)相当于把泰勒级数y(+~xx)=y()+(Ax)y+1(△)+…(—Ax)=y()一(△)y+1(△)+…截断于(Ax)v项,把(Ax)项以及更高幂次的项全部略去.式(3)相当于把泰勒级数y(+Ax)一Y(一△)=2(Ax)Y+(△)y,-+..截断于2(Ax)项,把(Ax)项以及更高幂次的项全部略去.因此,式(3)的误差小于式(1)和(2).二阶导数类似的可近似为差商的差商,一X[dx…一血dx【一止]:志[(+△)+y(—Ax)一2y()](4)[收稿日期]2011—11—21[作者简介]林喜季(1977一),女,福建江夏学院讲师,研究方向:代数表示论.16?第12卷?第1期林喜季:有限差分法在一维输运方程定解中的运用2012年2月这相当于把泰勒级数Y(+△)一v(一△)=2y()+(△)+(△)Y+..'截断于(△)项,把(△)项以及更高幂次的项全部略去.偏导数也可仿照式(1)一(4)近似为商差.这样一来,偏微分方程就成了差分方程.2一维输运方程的定解问题如,在区间(0,L)上求解一维输运方程"='axx.分析:(1)把整个空间分为.,个"步子",每一步的长度=I/J.于是,自变量以步长跳跃,它的取值是(i=0,1,2,…,.,).把时间步长取为zI,即自变量t取值t=kv(k=0,1,2,…,).(2)仿照式(1)和(4),一维输运方程可近似为(1)=(1—2lM(￡),2+l_!["(+l,t)+u(一1,t)】(5)这样只要知道某个时刻t的u在各个地点的值(,t),代人式(5)就可以得到下个时刻t…的的各个地点的值u(i,t).但这种解法时间t的步长z.不能太大,必须满足条件≤1,否则,由于舍入误差,会在其后各步的计算中产生雪崩影响,以致计算结果完全失去意义. (3)仿照式(2)和(4),一维输运方程可近似为u(,t)一U(i,t一1)r2(+1,t)+u(i一1,t)一2U(,t)即"(¨)=(1+2竿)Ⅱ(),2一旦j三[(+1,t)+"(,t)】(6)这样做可以取消对步长r的限制.但是知道某个时刻t的Ⅱ在各个地点的值(,t),并不能代入式(6)直接得到下个时刻t川的的各个地点的值(,t),且必须把i=1,2,3,…,.,一1的共计J一1个同式(6)的方程联立起来求解u(t,t+1),u(2,t+1),…,u(J一1,t%+1),当然这种联立方程的计算依靠电子计算机还是很方便的.(4)仿照式(3),偏导数近似为u(Xi'tk+1):,从而一维输运方程可近似为M(,t+1)一u(,t):u(Xi+?,tk+1)+u(—t,+.)一2u(,t+吉).上式中(,tk+)可理解为:[配(,t+.)+u(,t)】/2,于是,上例差分方程即为窘u(…)一(+字)"()+au("~i-1+)=一Ⅱ()一(1一窘)"(3Ci~tk)一骞H+1)(7)知道某个时刻t的在各个地点的值M(,t)后,必须把i=1,2,3,…,J一1的共计.,一1个同式(7)的方程联立起来求解"(.,t…),U(2,t+1),…,(J—l,t+1),当然这种联立方程的计算依靠电子计算机还是很方便的.这种解法对时间的步长也有限制,应满足2≤1,但与解法(2)相比限制要宽些.3用近似式求一维波动方程的定解问题如,在区间(0,L)上求解一维波动方程%一a"=0.把整个空间分为.,个"步子",每一步的长度=l/J.把时间步长取为r,仿照式(4),一维波动方程可近似为"(戈,t+1)+(,t一1)一2u(,t)即(.):2(1一譬)u(+旦}[(+1,tk)+(;一1,)一¨(,一1)】.这就是说,只要知道某时刻t及其以前时刻的U在各个地点i的值u(,t),代入上式,就可以得到下个时刻tk+.的的各个地点的值"(,t…).在此f青景下,时间的步长r的限制条件为≤1.17?,J,L2一,J一,L2"r,+,L22012年2月廊坊师范学院(自然科学版)第12卷?第1期从上述例子可以看出,在研究一维输运方程定解问题时,可采用有限差分法,按适当的数学变换把定解问题中的微商换成差商,从而把原问题离散化为差分格式,进而求出数值解.该方法具有简单,灵活,容易在计算机上实现的特点.并且该方法还有很强的通用性,如热传导过程,气体扩散过程这类定解问题,其过程都与时间有关,利用差分法解这类问题,就是从初始值出发,通过差分格式沿时间增加的方向,逐步求出微分方程的近似解.再如弹性力学中的平衡,电磁场及引力场等问题,其特征均为椭圆型方程,利用差分法解这类问题,就是合理选定的差分方格网,建立差分格式,最后求解代数方程组.[参考文献][1]吴顺唐,邓之光.有限差分法方程[M].南京:河海大学出版社.1993.[2]陈祖墀.偏微分方程[M].合肥:中国科学技术大学出版社.2004.[3]王晓东.算法与数据结构[M].北京:电子工业出版, 1998.[4]王震,谢树森.解四阶拟线性波动方程的一类二阶差分格式[J].中国海洋大学,2004,(34).(上接15页)将上述n为奇数与n为偶数两种情况统一起来,可得数列{a)的通项公式为.=1[(口.+.)+(一1)(.一n.)】?r.2)看P≠1,根据文献l3j结论司知,b=+()p,即有an+lan=+【一)p,从而anan_l=+(?一-qp)p.若令一1)=+(一,(n≥2),贝0ana一:f(一1).当p≠±,g≠o时,一)+(g一)p=一p棚=qp(1一p)≠0(n≥2),从而,('一p)≠0(n≥).由%a一.=/(n一1)递推可得:当n为大于1奇数时,(n一1)(一1)厂(//,一3).n':(—.—.::—0n一:je'.——.::—'':e';:—:二—;0n一=?…?f1al;一厂(n一2)厂(一4)厂(3)'()'当ll,为大于2的偶数时,n=:;{—;——{.一:=.一18?=船?…?n一12Ⅱ[,(2Jl})]即当n为大于1奇数时,a=丁_—一a.;Ⅱ厂()Ⅱ()当n为大于2的偶数时,.=T或改写为a=Ⅱ[,(2)]二者可统一为n=—1[(nl+口2)+(一1)(n2一口1)]?{{一吉【3+(一1)】)n一吉[3+(一1)】[Ⅱ[(2)]/Ⅱ厂(+)](-1),nEN+且n≥3,其中f(n)=(1一p),n∈N+.[参考文献][1]劳建祥.递推数列求通项大观[J].数学教学,2005,(3): 41—42.[2]高焕江.也谈二阶线性递推数列的周期?t:ff-[J].廊坊师范学院,2009,9(6):8—10.[3]高焕江.二阶线性递推数列的通项公式[J].保定学院学报,2010,23(3):34—37.。

Mixed Integer-Discrete-Continuous Optimization by Differential Evolution, Part 1 the optimi

Jouni Lampinen 1999
5th International Mendel Conference on Soft Computing MENDEL’99, Brno, Czech Republic, 9–12 June 1999
Jouni Lampinen – Ivan Zelinka Mixed Integer-Discrete-Continuous Optimization by Differential Evolution
trial vector cost value parameter 1 parameter 2 parameter 3 parameter 4 parameter 5 SELECTION: Select target vector or trial vector, the fittest one survive 2.56 -0.09 0.83 0.14 0.77 0.91
Байду номын сангаас
5th International Mendel Conference on Soft Computing MENDEL’99, Brno, Czech Republic, 9–12 June 1999
Jouni Lampinen 1999
Jouni Lampinen – Ivan Zelinka Mixed Integer-Discrete-Continuous Optimization by Differential Evolution
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 29 20 21 do 1 icount=1,GenMax do 2 i=1,NP A=IDINT(Rndgen()*NP+1) if(A.EQ.i)goto 3 B=IDINT(Rndgen()*NP+1) if((B.EQ.i).OR.(B.EQ.A))goto 4 C=IDINT(Rndgen()*NP+1) if((C.EQ.i).OR.(C.EQ.A).OR.(C.EQ.B))goto 5 do 6 j=1,iD trial(j)=P(i,j)+K*(P(C,j)-P(i,j))+F*(P(A,j)-P(B,j)) continue call CostFunction(score,iD,trial) if(score.LE.cost(i))then do 7 j=1,iD P(i,j)=trial(j) continue cost(i)=score iBest=i endif continue continue

abstract分析

Analysis of abstractAn abstract is a complete but concise statement that describes a larger work. Components of abstract are different from the discipline. The two most important reasons of abstracts are selection and indexing. Selection means that abstracts allow readers who may be interested in the long work to quickly decide whether this paper is worth their time to read it. Also, indexing means that many online databases use abstracts to index larger works. Therefore, the abstract is gaining tractions among us.There are two types of abstracts: Descriptive and Informative. A descriptive abstract outlines the topics covered in the paper. It makes no judgments about the work, nor does it provide results or conclusions of the research. The majority of abstracts are informative. An informative abstract includes the information that can be found in a descriptive abstract (purpose, methods, scope) but also includes the results and conclusions of the research and the recommendations of the author.From what I have mentioned above, the examples I choose (list of examples attached) can be divided into two: example 1 is a descriptive abstract, while the rests are all informative abstracts.The structure of a descriptive abstract may only include two parts like the example1: background and outline of the work. Usually, outline of the work likes a table of contents in paragraph form. Although descriptive abstracts are usually very short –100 words or less, it’s enough to know the structure of the whole work. Just li ke a miniature of the paper, it’s brief and complete.Example2~6 show the phenomenon that most abstracts are informative. An informative abstract often includes four parts: background, research question, method and conclusion. The background of an informative abstract is to answer the question that why do you want to do the research. Like example 2‖ Many practices aimed at cultivating multicultural competence in educational and organizational settings assume that multicultural experience fosters creativity‖, Is the assumption true? Or when and how does multicultural experience enhances creativity? With background of the topic, we will feel more comfortable to accept the conception of the work. The research question is the question the author trying to answer and the topic of the work. What problem you attempt to solve will determine whether this paper is worth the readers’ time.Like example 3‖This article investigates the possessions and activities that consumers love and their role in the construction of a coherent identity narrative‖, the question that example3 want to answer is the roles of loved objects and activities in the construction of a coherent identity narrative. The method of an informativeabstract is to answer the question that how do you get answers to your research question? While look at the list of examples, we’ll find that not all abstract show the approaches used in their work. The conclusion of an abstract is the work’s finding and value. Like example5 ―In an experiment, we show that vis ual personification—pictures in an ad that metaphorically represent a product as engaged in some kind of human behavior—can trigger anthropomorphism. Such personification, when embedded in an ad, appears to lead to more positive emotions, more positive attributions of brand personality, and greater brand liking. Implications for advertisers are discussed‖. From these words, we clearly know the final conclusion that visual personification can trigger anthropomorphism and personification in ads lead to positive effects.From the list of examples, we conclude that the language of abstract is brief and concise. Due to the limit of space, the abstract must be written the essential information in as simple terms as possible. Also, the work should be conveyed in as few words as possible with an assumption that it covers all the main point of the paper, that’s to say, the abstract should be complete. What’s more, the abstract ought to be reader-oriented, that’s means it should be written to serve reader in saving the ir time and energy.Here are some examples from the papers I have read recently:Example 1: Susan Fournier , Claudio Alvarez, Brands as relationship partners: Warmth, competence, and in-between, Journal of Consumer Psychology, October 2011.Background：The dialogue between social perception and consumer–brand relationship theories opens new opportunities for studying brands.Outline of the work:To advance branding research in the spirit of interdisciplinary inquiry, we propose to (1) investigate the process of anthropomorphism through which brands are imbued with intentional agency; (2) integrate the role of consumers not only as perceivers but also as relationship agents;(3) consider important defining dimensions of consumer–brand relationships beyond warmth and competence, including power and excitement; and (4) articulate the dynamics governing warmth (intentions) and competency (ability) judgments to yield prescriptive guidance for developing popular and admired brands.Example 2:Angela Ka-yee Leung, William W. Maddux ,Adam D. Galinsky, Chi-yue Chiu. Multicultural Experience Enhances Creativity:The When and How, American Psychologist, April 2008.Background：Many practices aimed at cultivating multicultural competence in educational and organizational settings (e.g., exchange programs, diversity education in college, diversity management at work) assume that multicultural experience fosters creativity.Research Question：In line with this assumption, the research reported in this article is the first to empirically demonstrate that exposure to multiple cultures in and of itself can enhance creativity.Conclusions:Overall, the authors found that extensiveness of multicultural experiences was positively related to both creative performance (insight learning, remote association, and idea generation) and creativity-supporting cognitive processes(retrieval of unconventional knowledge, recruitment of ideas from unfamiliar cultures for creative idea expansion). Furthermore, their studies showed that the serendipitous creative benefits resulting from multicultural experiences may depend on the extent to which individuals open themselves to foreign cultures, and that creativity is facilitated in contexts that deemphasize the need for firm answers or existential concerns.The authors discuss the implications of their findings for promoting creativity in increasingly global learning and work environments. Example 3: Aaron C. Ahuvia, Beyond the Extended Self: Loved Objects and Consumers’ Identity Narratives, Journal of Consumer Research, June 2005.Research Question：This article investigates the possessions and activities that consumers love and their role in the construction of a coherent identity narrative.Method：interviewConclusions: In the face of social forces pushing toward identity fragmentation, interviews reveal three different strategies, labeled ―demarcating,‖ ―compromising,‖ and ―synthesizing‖ solutions, for creating a coherent self-narrative.Findings are compared to Belk’s ―Possessions and the Extended Self.‖ Most claims from Belk are supported, but the notion of a core versus extended self is critiqued as a potentially confusing metaphor.The roles of loved objects and activities in structuring social relationships and in consumer well-being are also explored.Example 4:Nicolas Kervyn , Susan T. Fiske, Chris Malone, Brands as intentional agents framework: How perceived intentions and ability can map brand perception, Journal of Consumer Psychology, September 2011.Background：Building on the Stereotype Content Model, this paper introduces and tests the Brands as Intentional Agents Framework.A growing body of researchsuggests that consumers have relationships with brands that resemble relations between people.……Brands as Intentional Agents Framework is based on a well-established social perception approach: the Stereotype Content Model.Research Question：We propose that consumers perceive brands in the same way they perceive people. This approach allows us to explore how social perception theories and processes can predict brand purchase interest and loyalty.Conclusions: Two studies support the Brands as Intentional Agents Framework prediction that consumers assess a brand's perceived intentions and ability and that these perceptions elicit distinct emotions and drive differential brand behaviors.The research shows that human social interaction relationships translate to consumer–brand interactions in ways that are useful to inform brand positioning and brand communications.Example 5: Marjorie Delbaere, Edward F. McQuarrie, and Barbara J. Phillips, Personification in Advertising:Using a Visual Metaphor to Trigger Anthropomorphism, Journal of Advertising, Spring 2011.Background：All forms of personification draw on anthropomorphism, the propensity to attribute human characteristics to objects.Method：experiment,Conclusions: In an experiment, we show that visual personification—pictures in an ad that metaphorically represent a product as engaged in some kind of human behavior—can trigger anthropomorphism. Such personification, when embedded in an ad, appears to lead to more positive emotions, more positive attributions of brand personality, and greater brand liking. Implications for advertisers are discussed. Example 6:Marina Puzakova, Hyokjin Kwak,Joseph F. Rocereto, When Humanizing Brands Goes Wrong: The Detrimental Effect of Brand Anthropomorphization Amid Product Wrongdoings, Journal of Marketing, May 2013.Background：The brand relationship literature shows that the humanizing of brands and products generates more favorable consumer attitudes and thus enhances brand performance.Research Question：However, the authors propose negative downstream consequences of brand humanization; that is, the anthropomorphization of a brand can negatively affect consumers' brand evaluations when the brand faces negative publicity caused by product wrongdoings.Conclusions: They find that consumers who believe in personality stability (i.e., entity theorists) view anthropomorphized brands that undergo negative publicity less favorably than nonanthropomorphized brands. In contrast, consumers who advocate personality malleability (i.e.,incremental theorists) are less likely to devalue an anthropomorphized brand from a single instance of negative publicity.Finally, the authors explore three firm response strategies (i.e., denial,apology, and compensation) that can affect the evaluations of anthropomorphized brands for consumers with different implicit theory perspectives. They find that entity theorists have more difficulty in combating the adverse effects of brand anthropomorphization than incremental theorists. Furthermore, they demonstrate that compensation (vs. denial or apology) is the only effective response among entity theorists.。

Optimization Algorithms

Optimization AlgorithmsOptimization algorithms are a crucial tool in various fields, including engineering, economics, and computer science. These algorithms are designed tofind the best solution to a problem from a set of possible solutions, often with the goal of maximizing or minimizing a particular objective function. However, the effectiveness of optimization algorithms can vary depending on the specific problem at hand, and there are numerous factors to consider when selecting the most appropriate algorithm for a given situation. One of the key considerations when choosing an optimization algorithm is the nature of the problem itself. Some problems may be relatively simple and well-behaved, making them suitable for straightforward optimization techniques such as gradient descent or genetic algorithms. However, other problems may be highly complex, non-linear, or non-convex, requiring more advanced algorithms such as simulated annealing or particle swarm optimization. Understanding the characteristics of the problem is essential for selecting the most effective algorithm. Another important factor to consider is the computational resources available for running the optimization algorithm. Some algorithms may be highly computationally intensive, requiring significant time and memory resources to find a solution. In contrast, other algorithms may be more lightweight and suitable for use in resource-constrained environments. The availability of computational resources can significantly impact the choice of optimization algorithm. In addition to the nature of the problem and computational resources, the specific requirements and constraints of the problem must also be taken into account when selecting an optimization algorithm. For example, some problems may have strict constraints on the feasible solution space, requiring the use of constrained optimization algorithms. Likewise, certain problems may have multiple objectives that need to be simultaneously optimized, necessitating the use of multi-objective optimization algorithms. Understanding the specific requirements and constraints of the problem is crucial for choosing the most appropriate algorithm. Furthermore, the performance and robustness of optimization algorithms can vary depending on the specific problem and input parameters. Some algorithms may perform well on certain types of problems but poorly on others, while some algorithms may be sensitive to the choice ofparameters and initialization. It is essential to carefully evaluate the performance and robustness of different algorithms on a given problem to ensurethat the selected algorithm can reliably find high-quality solutions. Moreover, the interpretability and ease of implementation of optimization algorithms arealso important considerations. Some algorithms may produce highly complex and opaque solutions, making it challenging to understand the reasoning behind the optimized solution. In contrast, other algorithms may produce more interpretable solutions, which can be valuable in certain applications where explainability is critical. Additionally, the ease of implementation and integration of thealgorithm into existing systems and workflows should also be taken into account. Finally, the ethical and societal implications of using optimization algorithms should not be overlooked. In some cases, the use of optimization algorithms may have unintended consequences or ethical implications, such as reinforcing biasesor creating unfair outcomes. It is important to consider the potential ethical and societal impacts of using optimization algorithms and to take proactive measuresto mitigate any negative effects. In conclusion, the selection of an optimization algorithm is a complex and multifaceted decision that requires careful consideration of the problem at hand, computational resources, specific requirements and constraints, algorithm performance and robustness,interpretability and ease of implementation, and ethical and societal implications. By taking into account these various factors, practitioners can make informed decisions about which optimization algorithm is best suited for a given problem, ultimately leading to more effective and responsible use of optimization techniques.。

abstract单词讲解

abstract单词讲解"Abstract" is an adjective that means existing in thought or as an idea but not having a physical or concrete existence. It can also be used as a verb, meaning to consider something in a general or non-specific way.In the context of art or design, "abstract" refers to a style that does not attempt to represent an accurate depiction of visual reality. Instead, it focuses on shapes, colors, forms, and gestural marks to achieve its effect. Abstract art often aims to evoke emotions, sensations, or ideas, rather than portraying recognizable objects or scenes.In academic writing, an "abstract" is a brief summary of a research paper, typically found at the beginning of the paper. It provides an overview of the study's objectives, methods, results, and conclusions. The purpose of an abstract is to allow readers to quickly determine whether or not they are interested in reading the full paper.In general, the term "abstract" is used to describe something that is conceptual, theoretical, or non-physical in nature. It is often contrasted with "concrete" or "literal" concepts, which refer to things that can be directly perceived or experienced.。

复杂产品概念设计多学科过程建模方法研究

刘海强，祁国宁，张太华,纪杨建（浙江大学现代制造工程研究所email:liuhaiqiang1980@）摘要：为了使复杂产品概念设计的过程模型支持多学科优化，实现概念设计过程中的求解功能，确定产品设计的综合最优方案，提出了基于扩展Petri网与NSGA-Ⅱ算法相结合的多学科过程建模方法。

针对复杂产品设计过程的特点，分析了概念设计在多学科优化过程中的重要作用，研究了扩展Petri网模型建立和NSGA-Ⅱ算法的执行过程。

重点讨论以扩展Petri模型的变迁序列作为NSGA-Ⅱ算法的染色体,结合属于Petri网模型元素的选择算子、交叉算子进行NSGA-Ⅱ优化运算。

将上述方法应用于某工业汽轮机的设计，可以对于不同的设计要求，给出不同情况下的Pareto最优解集，是进行复杂产品概念设计过程建模的一种切实有效方法。

关键词：过程建模；概念设计；多学科；Pareto最优中图分类号：TH122复杂产品是指结构复杂、技术密集、制造要求和成本高、过程管理复杂、客户需求多变的一类产品[1]。

复杂产品的设计是一个跨越多学科的、复杂且高度迭代的过程，这个过程涉及到多专业信息的集成、多学科团队的协同设计，其复杂性主要表现在各种设计因素（如设计零部件、设计活动等）之间的相互耦合作用。

衡量一个产品的设计是否成功标准有不少，取决于用户所最追求的目标，这个目标通常是产品成本、周期、质量等。

概念设计所花费的成本只占整个产品全生命周期成本的约1％，但它却决定了整个全生命周期成本的70％[2]。

概念设计不仅决定着开发产品的成本、周期、质量、性能、可靠性、安全性和环保性，而且产生的设计无法由后续设计过程弥补。

概念设计被认为是设计过程中最重要、最关键、最具创造性的阶段[3]。

概念设计对产品全生命周期的其它环节，如生产、制造、维修等也有着重要的影响。

过程建模是对产品开发过程的表示和设计,它是产品开发过程研究和应用的首要问题。

过程模型是产品开发过程及其内在规律的一种形式化描述，是对开发设计过程本质的一种表示。

Ant algorithms for discrete optimization

Ant Algorithms for Discrete Optimization Marco DorigoGianni Di CaroIRIDIA CP194/6Universit´e Libre de Bruxelles Avenue Franklin Roosevelt50 B-1050BrusselsBelgiummdorigo@ulb.ac.begdicaro@iridia.ulb.ac.beLuca M.GambardellaIDSIACorso Elvezia36CH-6900LuganoSwitzerlandluca@idsia.chKeywordsant algorithms,ant colony optimiza-tion,swarm intelligence,metaheuris-tics,natural computationAbstract This article presents an overview of recent workon ant algorithms,that is,algorithms for discrete optimizationthat took inspiration from the observation of ant colonies’foraging behavior,and introduces the ant colony optimization(ACO)metaheuristic.In theﬁrst part of the article the basicbiologicalﬁndings on real ants are reviewed and theirartiﬁcial counterparts as well as the ACO metaheuristic aredeﬁned.In the second part of the article a number ofapplications of ACO algorithms to combinatorial optimizationand routing in communications networks are described.Weconclude with a discussion of related work and of some ofthe most important aspects of the ACO metaheuristic.1IntroductionAnt algorithms wereﬁrst proposed by Dorigo and colleagues[33,40]as a multi-agent approach to difﬁcult combinatorial optimization problems such as the traveling sales-man problem(TSP)and the quadratic assignment problem(QAP).There is currently much ongoing activity in the scientiﬁc community to extend and apply ant-based al-gorithms to many different discrete optimization problems[5,21].Recent applications cover problems such as vehicle routing,sequential ordering,graph coloring,routing in communications networks,and so on.Ant algorithms were inspired by the observation of real ant colonies.Ants are social insects,that is,insects that live in colonies and whose behavior is directed more to the survival of the colony as a whole than to that of a single individual component of the colony.Social insects have captured the attention of many scientists because of the high structuration level their colonies can achieve,especially when compared to the relative simplicity of the colony’s individuals.An important and interesting behavior of ant colonies is their foraging behavior,and,in particular,how ants canﬁnd the shortest paths between food sources and their nest.While walking from food sources to the nest and vice versa,ants deposit on the ground a substance called pheromone,forming in this way a pheromone trail.Ants can smell pheromone,and when choosing their way,they tend to choose,in probability, paths marked by strong pheromone concentrations.The pheromone trail allows the ants toﬁnd their way back to the food source(or to the nest).Also,it can be used by other ants toﬁnd the location of the food sources found by their nestmates.It has been shown experimentally that this pheromone trail following behavior can give rise,once employed by a colony of ants,to the emergence of shortest paths.That is,when more paths are available from the nest to a food source,a colony of ants may be able to exploit the pheromone trails left by the individual ants to discover the shortest path from the nest to the food source and back.c 1999Massachusetts Institute of Technology Artiﬁcial Life5:137–172(1999)(a)Time (minutes)(a)(b)020*********1The above-described experiments have been run in strongly constrained conditions.A formal proof of the pheromone-driven shortest path-ﬁnding behavior in the general case is missing.Bruckstein et al.[9,10]consider the shortest path-ﬁnding problem in absence of obstacles for ants driven by visual clues and not by pheromones and prove the convergence of the ants’path to the straight line.Artiﬁcial Life Volume5,Number2139As deﬁned by Grass´e in his work on Bellicositermes Natalensis and Cubitermes [62],stigmergy is the“stimulation of the workers by the very performances they have achieved”(p.79).2In fact,Grass´e[61]observed that insects are capable of responding to so-called “signiﬁcant stimuli”that activate a genetically encoded reaction.In social insects,of which termites and ants are some of the best known examples,the effects of these reactions can act as new signiﬁcant stimuli for both the insect that produced them and for other insects in the colony.The production of a new signiﬁcant stimulus as a consequence of the reaction to a signiﬁcant stimulus determines a form of coordination of the activities and can be interpreted as a form of indirect communication.For example,Grass´e[62]observed that Bellicositermes Natalensis as well as Cubitermes, when building a new nest,start by a random,noncoordinated activity of earth pellet depositing.But once the earth pellets reach a certain density in a restricted area they become a new signiﬁcant stimulus that causes more termites to add earth pellets so that pillar and arches,and eventually the whole nest,are built.What characterizes stigmergy from other means of communication is(a)the physical nature of the information released by the communicating insects,which corresponds to a modiﬁcation of physical environmental states visited by the insects,and(b)the local nature of the released information,which can only be accessed by insects that visit the state in which it was released(or some neighborhood of that state).Accordingly,in this article we take the stance that it is possible to talk of stig-mergetic communication whenever there is an indirect communication mediated by physical modiﬁcations of environmental states which are only locally accessible by the communicating agents.One of the main tenets of this article is that the stigmergetic model of communi-cation in general,and the one inspired by ants’foraging behavior in particular,is an interesting model for artiﬁcial multi-agent systems applied to the solution of difﬁcult optimization problems.In fact,the above-mentioned characteristics of stigmergy can easily be extended to artiﬁcial agents by(a)associating with problem states appropriate state variables,and(b)giving the artiﬁcial agents only local access to these variables’values.For example,in the above-described foraging behavior of ants,stigmergetic commu-nication is at work via the pheromone that ants deposit on the ground while walking. Correspondingly,our artiﬁcial ants will simulate pheromone laying by modifying appro-priate“pheromone variables”associated with problem states they visit while building solutions to the optimization problem to which they are applied.Also,according to the stigmergetic communication model,our artiﬁcial ants will have only local access to these pheromone variables.Another important aspect of real ants’foraging behavior that is exploited by artiﬁcial ants is the coupling between the autocatalytic(positive feedback)mechanism[40]and the implicit evaluation of solutions.By implicit solution evaluation we mean the fact that shorter paths(which correspond to lower cost solutions in artiﬁcial ants)will be completed earlier than longer ones,and therefore they will receive pheromone rein-forcement more quickly.Implicit solution evaluation coupled with autocatalysis can be very effective:The shorter the path,the sooner the pheromone is deposited by the ants,the more the ants that use the shorter path.If appropriately used,autocatalysis can be a powerful mechanism in population-based optimization algorithms(e.g.,in evolutionary computation algorithms[45,66,85,91]autocatalysis is implemented by the selection/reproduction mechanism).In fact,it quickly favors the best individuals, 2Workers are one of the castes in termite colonies.Although Grass´e introduced the term stigmergy to explain the behavior of termite societies,the same term has been used to describe indirect communication mediated by modiﬁcations of the environment that can also be observed in other social insects.140Artiﬁcial Life Volume5,Number2so that they can direct the search process.When using autocatalysis some care must be taken to avoid premature convergence(stagnation),that is,the situation in which some not very good individual takes over the population just because of a contingent situation(e.g.,because of a local optimum,or just because of initial randomﬂuctu-ations that caused a not very good individual to be much better than all the other individuals in the population)impeding further exploration of the search space.We will see that pheromone trail evaporation and stochastic state transitions are the needed complements to autocatalysis drawbacks.In the remainder of this article we discuss a number of ant algorithms based on the above ideas.We start by deﬁning,in Section2,the characterizing aspects of ant algorithms and the ant colony optimization(ACO)metaheuristic.3Section3is an overview of most of the applications of ACO algorithms.In Section4we brieﬂy discuss related work,and in Section5we discuss some of the characteristics of implemented ACO algorithms.Finally,we draw some conclusions in Section6.2The Ant Colony Optimization ApproachIn the ant colony optimization(ACO)metaheuristic a colony of artiﬁcial ants cooperates inﬁnding good solutions to difﬁcult discrete optimization problems.Cooperation is a key design component of ACO algorithms:The choice is to allocate the computational resources to a set of relatively simple agents(artiﬁcial ants)that communicate indirectly by stigmergy.Good solutions are an emergent property of the agents’cooperative interaction.Artiﬁcial ants have a double nature.On the one hand,they are an abstraction of those behavioral traits of real ants that seemed to be at the heart of the shortest path-ﬁnding behavior observed in real ant colonies.On the other hand,they have been enriched with some capabilities that do notﬁnd a natural counterpart.In fact,we want ant colony optimization to be an engineering approach to the design and implementation of software systems for the solution of difﬁcult optimization problems.It is therefore reasonable to give artiﬁcial ants some capabilities that,although not corresponding to any capacity of their real ant counterparts,make them more effective and efﬁcient.In the following we discussﬁrst the nature-inspired characteristics of artiﬁcial ants,and then how they differ from real ants.2.1Similarities and Differences with Real AntsMost of the ideas of ACO stem from real ants.In particular,the use of:(a)a colony of cooperating individuals,(b)an(artiﬁcial)pheromone trail for local stigmergetic com-munication,(c)a sequence of local moves toﬁnd shortest paths,and(d)a stochastic decision policy using local information and no lookahead.Colony of cooperating individuals.As real ant colonies,ant algorithms are composed of a population,or colony,of concurrent and asynchronous entities globally cooperat-ing toﬁnd a good“solution”to the task under consideration.Although the complexity of each artiﬁcial ant is such that it can build a feasible solution(as a real ant can some-howﬁnd a path between the nest and the food),high quality solutions are the result of the cooperation among the individuals of the whole colony.Ants cooperate by means of the information they concurrently read/write on the problem’s states they visit,as explained in the next item.3It is important here to clarify brieﬂy the terminology used.We talk of ACO metaheuristic to refer to the general procedure presented in Section2.The term ACO algorithm will be used to indicate any generic instantiation of the ACO metaheuristic. Alternatively,we will also talk more informally of ant algorithms to indicate any algorithm that,while following the general guidelines set above,does not necessarily follow all the aspects of the ACO metaheuristic.Therefore,all ACO algorithms are also ant algorithms,though the converse is not true(e.g.,we will see that HAS-QAP is an ant,but not strictly an ACO algorithm). Artiﬁcial Life Volume5,Number2141Pheromone trail and stigmergy.Artiﬁcial ants modify some aspects of their envi-ronment as the real ants do.While real ants deposit on the world’s state they visit a chemical substance,the pheromone,artiﬁcial ants change some numeric information locally stored in the problem’s state they visit.This information takes into account the ant’s current history or performance and can be read/written by any ant access-ing the state.By analogy,we call this numeric information artiﬁcial pheromone trail, pheromone trail for short in the following.In ACO algorithms local pheromone trails are the only communication channels among the ants.This stigmergetic form of com-munication plays a major role in the utilization of collective knowledge.Its main effect is to change the way the environment(the problem landscape)is locally perceived by the ants as a function of all the past history of the whole ant ually,in ACO algorithms an evaporation mechanism,similar to real pheromone evaporation,modi-ﬁes pheromone information over time.Pheromone evaporation allows the ant colony slowly to forget its past history so that it can direct its search toward new directions without being over-constrained by past decisions.Shortest path searching and local moves.Artiﬁcial and real ants share a common task:toﬁnd a shortest(minimum cost)path joining an origin(nest)to destination (food)sites.Real ants do not jump;they just walk through adjacent terrain’s states,and so do artiﬁcial ants,moving step-by-step through“adjacent states”of the problem.Of course,exact deﬁnitions of state and adjacency are problem speciﬁc.Stochastic and myopic state transition policy.Artiﬁcial ants,as real ones,build so-lutions applying a probabilistic decision policy to move through adjacent states.As for real ants,the artiﬁcial ants’policy makes use of local information only and it does not make use of lookahead to predict future states.Therefore,the applied policy is com-pletely local,in space and time.The policy is a function of both the a priori information represented by the problem speciﬁcations(equivalent to the terrain’s structure for real ants),and of the local modiﬁcations in the environment(pheromone trails)induced by past ants.As we said,artiﬁcial ants also have some characteristics that do notﬁnd their coun-terpart in real ants.•Artiﬁcial ants live in a discrete world and their moves consist of transitions from discrete states to discrete states.•Artiﬁcial ants have an internal state.This private state contains the memory of the ants’past actions.•Artiﬁcial ants deposit an amount of pheromone that is a function of the quality of the solution found.4•Artiﬁcial ants’timing in pheromone laying is problem dependent and often does not reﬂect real ants’behavior.For example,in many cases artiﬁcial ants update pheromone trails only after having generated a solution.•To improve overall system efﬁciency,ACO algorithms can be enriched with extra capabilities such as lookahead,local optimization,backtracking,and so on that cannot be found in real ants.In many implementations ants have been hybridized with local optimization procedures(see,e.g.,[38,51,98]),while,so far,only Michel and Middendorf[78]have used a simple one-step lookahead function and there are no examples of backtracking procedures added to the basic ant capabilities,except for simple recovery procedures used by Di Caro and Dorigo[26,29].54In reality,some real ants have a similar behavior:They deposit more pheromone in case of richer food sources.5Usually,backtracking strategies are suitable to solve constraint satisfaction problems,(e.g.,n-queens)and lookahead is very useful 142Artiﬁcial Life Volume5,Number2In the following section we will show how artiﬁcial ants can be put to work in an algorithmic framework so that they can be applied to discrete optimization problems.2.2The ACO MetaheuristicIn ACO algorithms aﬁnite-size colony of artiﬁcial ants with the above-described char-acteristics collectively searches for good-quality solutions to the optimization prob-lem under consideration.Each ant builds a solution,or a component of it,6starting from an initial state selected according to some problem-dependent criteria.While building its own solution,each ant collects information on the problem character-istics and on its own performance and uses this information to modify the repre-sentation of the problem,as seen by the other ants.Ants can act concurrently and independently,showing a cooperative behavior.They do not use direct communica-tion:It is the stigmergy paradigm that governs the information exchange among the ants.An incremental constructive approach is used by the ants to search for a feasi-ble solution.A solution is expressed as a minimum cost(shortest)path through the states of the problem in accordance with the problem’s constraints.The com-plexity of each ant is such that even a single ant is able toﬁnd a(probably poor quality)solution.High-quality solutions are only found as the emergent result of the global cooperation among all the agents of the colony concurrently building different solutions.According to the assigned notion of neighborhood(problem-dependent),each ant builds a solution by moving through a(ﬁnite)sequence of neighbor states.Moves are selected by applying a stochastic local search policy directed(a)by ant private information(the ant internal state,or memory)and(b)by publicly available pheromone trails and a priori problem-speciﬁc local information.The ant’s internal state stores information about the ant past history.It can be used to carry useful information to compute the value/goodness of the generated so-lution and/or the contribution of each executed move.Moreover it can play a fun-damental role to manage the feasibility of the solutions.In some problems,in fact, typically in combinatorial optimization,some of the moves available to an ant in a state can take the ant to an infeasible state.This can be avoided by exploiting the ant’s memory.Ants therefore can build feasible solutions using only knowledge about the local state and about the effects of actions that can be performed in the local state.The local,public information comprises both some problem-speciﬁc heuristic infor-mation and the knowledge,coded in the pheromone trails,accumulated by all the ants from the beginning of the search process.This time-global pheromone knowledge built up by the ants is a shared local long-term memory that inﬂuences the ants’decisions. The decisions about when the ants should release pheromone on the“environment”and how much pheromone should be deposited depend on the characteristics of the problem and on the design of the implementation.Ants can release pheromone while building the solution(online step-by-step),or after a solution has been built,moving back to all the visited states(online delayed),or both.As we said,autocatalysis plays when the cost of making a local prediction about the effect of future moves is much lower than the cost of the real execution of the move sequence(e.g.,mobile robotics).T o our knowledge,until now ACO algorithms have not been applied to these classes of problems.6T o make more intuitive what we mean by component of a solution,we can consider,as an example,a transportation routing problem:Given a set of n cities,{c i},i∈{1,...,n},and a network of interconnection roads,we want toﬁnd all the shortest paths s ij connecting each city pair c i c j.In this case,a complete solution is represented by the set of all the n(n−1)shortest path pairs, while a component of a solution is a single path s ij.Artiﬁcial Life Volume5,Number2143an important role in ACO algorithms functioning:The more ants choose a move,the more the move is rewarded(by adding pheromone)and the more interesting it be-comes for the next ants.In general,the amount of pheromone deposited is made proportional to the goodness of the solution an ant has built(or is building).In this way,if a move contributed to generating a high-quality solution its goodness will be increased proportionally to its contribution.A functional composition of the locally available pheromone and heuristic values deﬁnes ant-decision tables,that is,probabilistic tables used by the ants’decision policy to direct their search toward the most interesting regions of the search space.The stochastic component of the move choice decision policy and the previously dis-cussed pheromone evaporation mechanism prevent a rapid drift of all the ants to-ward the same part of the search space.Of course,the level of stochasticity in the policy and the strength of the updates in the pheromone trail determine the bal-ance between the exploration of new points in the state space and the exploitation of accumulated knowledge.If necessary and feasible,the ants’decision policy can be enriched with problem-speciﬁc components such as backtracking procedures or lookahead.Once an ant has accomplished its task,consisting of building a solution and depositing pheromone information,the ant“dies,”that is,it is deleted from the system.The overall ACO metaheuristic,besides the two above-described components acting from a local perspective(i.e.,ants’generation and activity,and pheromone evapora-tion),can also comprise some extra components that use global information and that go under the name of daemon actions in the algorithm reported in Figure3.For example,a daemon can be allowed to observe the ants’behavior and to collect useful global information to deposit additional pheromone information,biasing,in this way, the ant search process from a nonlocal perspective.Or,it could,on the basis of the observation of all the solutions generated by the ants,apply problem-speciﬁc local op-timization procedures and deposit additional pheromone“ofﬂine”with respect to the pheromone the ants deposited online.The three main activities of an ACO algorithm(ant generation and activity,pheromone evaporation,and daemon actions)may need some kind of synchronization,performed by the schedule activities construct of Figure3.In general,a strictly sequential scheduling of the activities is particularly suitable for nondistributed problems,where the global knowledge is easily accessible at any instant and the operations can be conveniently synchronized.On the contrary,some form of parallelism can be easily and efﬁciently exploited in distributed problems such as routing in telecommunications networks,as will be discussed in Section3.2.In Figure3,a high-level description of the ACO metaheuristic is reported in pseudo-code.As pointed out above,some described components and behaviors are optional, such as daemon activities,or strictly implementation dependent,such as when and how the pheromone is deposited.In general,the online step-by-step pheromone update and the online delayed pheromone update components(respectively,lines24–27and 30–34in the new active ant()procedure)are mutually exclusive and only in a few cases are they both present or both absent(when both components are absent,the pheromone is deposited by the daemon).ACO algorithms,as a consequence of their concurrent and adaptive nature,are particularly suitable for distributed stochastic problems where the presence of exoge-nous sources determines a nonstationarity in the problem representation(in terms of costs and/or environment).For example,many problems related to communications or transportation networks are intrinsically distributed and nonstationary and it is often not possible to have an exact model of the underlying variability.On the contrary, because stigmergy is both the only inter-ant communication method and it is spatially 144Artiﬁcial Life Volume5,Number2Figure3.The ACO metaheuristic in ments are enclosed in braces.All the procedures at theﬁrst level of indentation in the statement in parallel are executed concurrently.The procedure daemon actions() at line6is optional and refers to centralized actions executed by a daemon possessing global knowledge.The target state(line19)refers to a complete solution,or to a component of a complete solution,built by the ant.The step-by-step and delayed pheromone updating procedures at lines24–27and30–34are often mutually exclusive.When both of them are absent the pheromone is deposited by the daemon.localized,ACO algorithms could perform not at their best in problems where each state has a large-sized neighborhood.In fact,an ant that visits a state with a large-sized neighborhood has a huge number of possible moves among which to choose.There-fore,the probability that many ants visit the same state is very small,and consequently there is little,if any,difference between using or not using pheromone trails.A more formal deﬁnition of the ACO metaheuristic,as well as of the class of problems to which it can be applied,can be found in[35].Artiﬁcial Life Volume5,Number21453Applications of ACO AlgorithmsThere are now available numerous successful implementations of the ACO metaheuris-tic(Figure3)applied to a number of different combinatorial optimization problems. Looking at these implementations it is possible to distinguish among two classes of ap-plications:those to static combinatorial optimization problems,and those to dynamic ones.Static problems are those in which the characteristics of the problem are given once and for all when the problem is deﬁned and do not change while the problem is being solved.A paradigmatic example of such problems is the classic traveling salesman problem[67,71,86],in which city locations and their relative distances are part of the problem deﬁnition and do not change at run time.On the contrary,dynamic problems are deﬁned as a function of some quantities whose value is set by the dynamics of an underlying system.The problem changes therefore at run time and the optimiza-tion algorithm must be capable of adapting online to the changing environment.The paradigmatic example discussed in the remainder of this section is network routing.Topological modiﬁcations(e.g.,adding or removing a node),which are not consid-ered by the above classiﬁcation,can be seen as transitions between problems belonging to the same class.Tables1and2list the available implementations of ACO algorithms.The main characteristics of the listed algorithms are discussed in the following two subsections. We then conclude with a brief review of existing parallel implementations of ACO algorithms.3.1Applications of ACO Algorithms to Static Combinatorial OptimizationProblemsThe application of the ACO metaheuristic to a static combinatorial optimization problem is relatively straightforward,once a mapping of the problem that allows the incremental construction of a solution,a neighborhood structure,and a stochastic state transition rule to be locally used to direct the constructive procedure is deﬁned.A strictly implementation-dependent aspect of the ACO metaheuristic regards the timing of pheromone updates(lines24–27and30–34of the algorithm in Figure3).In ACO algorithms for static combinatorial optimization the way ants update pheromone trails changes across algorithms:Any combination of online step-by-step pheromone updates and online delayed pheromone updates is possible.Another important implementation-dependent aspect concerns the daemon actions()component of the ACO metaheuristic(line6of the algorithm in Figure3). Daemon actions implement actions that require some kind of global knowledge about the problem.Examples are ofﬂine pheromone updates and local optimization of solu-tions built by ants.Most of the ACO algorithms presented in this subsection are strongly inspired by Ant System(AS),theﬁrst work on ant colony optimization[33,40].Many of the successive applications of the original idea are relatively straightforward applications of AS to the speciﬁc problem under consideration.We therefore start the description of ACO algorithms with AS.Following AS,for each ACO algorithm listed in Table1 we give a short description of the algorithm’s main characteristics and of the results obtained.3.1.1T raveling Salesman ProblemTheﬁrst application of an ant colony optimization algorithm was done using the trav-eling salesman problem as a test problem.The main reasons why the TSP,one of the most studied NP-hard[71,86]problems in combinatorial optimization,was chosen are 146Artiﬁcial Life Volume5,Number2。

Abstract

DMM’s topology is detected automatically dependent on execution of a special program. This paper outlines the mpC features concerning networks and presents the implementation of network management in the mpC programming environment. Details of the language are presented elsewhere [14mpC language was developed to write efﬁcient and portable programs for wide range of distributed memory machines. It supports both task and data parallelism, allows both static and dynamic process and communication structures, enables optimizations aimed at both communication and computation, and supports modular parallel programming and the development of a library of parallel programs. The language is an ANSI C superset based on the notion of network comprising processor nodes of different types and performances, connected with links of different bandwidths. The user can describe a network topology, create and discard networks, distribute data and computations over networks. The mpC programming environment uses the topological information in run time to ensure the efﬁcient execution of the application.The paper describes the implementation of managing networks in the mpC programming environment.

Multi-objective Optimization

Chapter2Multi-objective OptimizationAbstract In this chapter,we introduce multi-objective optimization,and recall some of the most relevant research articles that have appeared in the international litera-ture related to these topics.The presented state-of-the-art does not have the purpose of being exhaustive;it aims to drive the reader to the main problems and the ap-proaches to solve them.2.1Multi-objective ManagementThe choice of a route at a planning level can be done taking into account time, length,but also parking or maintenance facilities.As far as advisory or,more in general,automation procedures to support this choice are concerned,the available tools are basically based on the“shortest-path problem”.Indeed,the problem toﬁnd the single-objective shortest path from an origin to a destination in a network is one of the most classical optimization problems in transportation and logistic,and has deserved a great deal of attention from researchers worldwide.However,the need to face real applications renders the hypothesis of a single-objective function to be optimized subject to a set of constraints no longer suitable,and the introduction of a multi-objective optimization framework allows one to manage more informa-tion.Indeed,if for instance we consider the problem to route hazardous materials in a road network(see,e.g.,Erkut et al.,2007),deﬁning a single-objective function problem will involve,separately,the distance,the risk for the population,and the transportation costs.If we regard the problem from different points of view,i.e.,in terms of social needs for a safe transshipment,or in terms of economic issues or pol-11122Multi-objective Optimizationlution reduction,it is clear that a model that considers simultaneously two or more such objectives could produce solutions with a higher level of equity.In the follow-ing,we will discuss multi-objective optimization and related solution techniques.2.2Multi-objective Optimization and Pareto-optimal SolutionsA basic single-objective optimization problem can be formulated as follows:min f(x)x∈S,where f is a scalar function and S is the(implicit)set of constraints that can be deﬁned asS={x∈R m:h(x)=0,g(x)≥0}.Multi-objective optimization can be described in mathematical terms as follows:min[f1(x),f2(x),...,f n(x)]x∈S,where n>1and S is the set of constraints deﬁned above.The space in which the objective vector belongs is called the objective space,and the image of the feasible set under F is called the attained set.Such a set will be denoted in the following withC={y∈R n:y=f(x),x∈S}.The scalar concept of“optimality”does not apply directly in the multi-objective setting.Here the notion of Pareto optimality has to be introduced.Essentially,a vector x∗∈S is said to be Pareto optimal for a multi-objective problem if all other vectors x∈S have a higher value for at least one of the objective functions f i,with i=1,...,n,or have the same value for all the objective functions.Formally speak-ing,we have the following deﬁnitions:•A point x∗is said to be a weak Pareto optimum or a weak efﬁcient solution for the multi-objective problem if and only if there is no x∈S such that f i(x)<f i(x∗) for all i∈{1,...,n}.2.2Multi-objective Optimization and Pareto-optimal Solutions13•A point x∗is said to be a strict Pareto optimum or a strict efﬁcient solution for the multi-objective problem if and only if there is no x∈S such that f i(x)≤f i(x∗) for all i∈{1,...,n},with at least one strict inequality.We can also speak of locally Pareto-optimal points,for which the deﬁnition is the same as above,except that we restrict attention to a feasible neighborhood of x∗.In other words,if B(x∗,ε)is a ball of radiusε>0around point x∗,we require that for someε>0,there is no x∈S∩B(x∗,ε)such that f i(x)≤f i(x∗)for all i∈{1,...,n}, with at least one strict inequality.The image of the efﬁcient set,i.e.,the image of all the efﬁcient solutions,is called Pareto front or Pareto curve or surface.The shape of the Pareto surface indicates the nature of the trade-off between the different objective functions.An example of a Pareto curve is reported in Fig.2.1,where all the points between(f2(ˆx),f1(ˆx))and (f2(˜x),f1(˜x))deﬁne the Pareto front.These points are called non-inferior or non-dominated points.f1(xFig.2.1Example of a Pareto curveAn example of weak and strict Pareto optima is shown in Fig.2.2:points p1and p5are weak Pareto optima;points p2,p3and p4are strict Pareto optima.142Multi-objective Optimization2Fig.2.2Example of weak and strict Pareto optima2.3Techniques to Solve Multi-objective Optimization ProblemsPareto curves cannot be computed efﬁciently in many cases.Even if it is theoreti-cally possible toﬁnd all these points exactly,they are often of exponential size;a straightforward reduction from the knapsack problem shows that they are NP-hard to compute.Thus,approximation methods for them are frequently used.However, approximation does not represent a secondary choice for the decision maker.Indeed, there are many real-life problems for which it is quite hard for the decision maker to have all the information to correctly and/or completely formulate them;the deci-sion maker tends to learn more as soon as some preliminary solutions are available. Therefore,in such situations,having some approximated solutions can help,on the one hand,to see if an exact method is really required,and,on the other hand,to exploit such a solution to improve the problem formulation(Ruzica and Wiecek, 2005).Approximating methods can have different goals:representing the solution set when the latter is numerically available(for convex multi-objective problems);ap-proximating the solution set when some but not all the Pareto curve is numerically available(see non-linear multi-objective problems);approximating the solution set2.3Techniques to Solve Multi-objective Optimization Problems15when the whole efﬁcient set is not numerically available(for discrete multi-objective problems).A comprehensive survey of the methods presented in the literature in the last33 years,from1975,is that of Ruzica and Wiecek(2005).The survey analyzes sepa-rately the cases of two objective functions,and the case with a number of objective functions strictly greater than two.More than50references on the topic have been reported.Another interesting survey on these techniques related to multiple objec-tive integer programming can be found in the book of Ehrgott(2005)and the paper of Erghott(2006),where he discusses different scalarization techniques.We will give details of the latter survey later in this chapter,when we move to integer lin-ear programming formulations.Also,T’Kindt and Billaut(2005)in their book on “Multicriteria scheduling”,dedicated a part of their manuscript(Chap.3)to multi-objective optimization approaches.In the following,we will start revising,following the same lines of Erghott (2006),these scalarization techniques for general continuous multi-objective op-timization problems.2.3.1The Scalarization TechniqueA multi-objective problem is often solved by combining its multiple objectives into one single-objective scalar function.This approach is in general known as the weighted-sum or scalarization method.In more detail,the weighted-sum method minimizes a positively weighted convex sum of the objectives,that is,minn∑i=1γi·f i(x)n∑i=1γi=1γi>0,i=1,...,nx∈S,that represents a new optimization problem with a unique objective function.We denote the above minimization problem with P s(γ).It can be proved that the minimizer of this single-objective function P(γ)is an efﬁcient solution for the original multi-objective problem,i.e.,its image belongs to162Multi-objective Optimizationthe Pareto curve.In particular,we can say that if theγweight vector is strictly greater than zero(as reported in P(γ)),then the minimizer is a strict Pareto optimum,while in the case of at least oneγi=0,i.e.,minn∑i=1γi·f i(x)n∑i=1γi=1γi≥0,i=1,...,nx∈S,it is a weak Pareto optimum.Let us denote the latter problem with P w(γ).There is not an a-priori correspondence between a weight vector and a solution vector;it is up to the decision maker to choose appropriate weights,noting that weighting coefﬁcients do not necessarily correspond directly to the relative impor-tance of the objective functions.Furthermore,as we noted before,besides the fact that the decision maker cannot be aware of which weights are the most appropriate to retrieve a satisfactorily solution,he/she does not know in general how to change weights to consistently change the solution.This means also that it is not easy to develop heuristic algorithms that,starting from certain weights,are able to deﬁne iteratively weight vectors to reach a certain portion of the Pareto curve.Since setting a weight vector conducts to only one point on the Pareto curve,per-forming several optimizations with different weight values can produce a consid-erable computational burden;therefore,the decision maker needs to choose which different weight combinations have to be considered to reproduce a representative part of the Pareto front.Besides this possibly huge computation time,the scalarization method has two technical shortcomings,as explained in the following.•The relationship between the objective function weights and the Pareto curve is such that a uniform spread of weight parameters,in general,does not producea uniform spread of points on the Pareto curve.What can be observed aboutthis fact is that all the points are grouped in certain parts of the Pareto front, while some(possibly signiﬁcative)portions of the trade-off curve have not been produced.2.3Techniques to Solve Multi-objective Optimization Problems17•Non-convex parts of the Pareto set cannot be reached by minimizing convex combinations of the objective functions.An example can be made showing a geometrical interpretation of the weighted-sum method in two dimensions,i.e., when n=2.In the two-dimensional space the objective function is a liney=γ1·f1(x)+γ2·f2(x),wheref2(x)=−γ1·f1(x)γ2+yγ2.The minimization ofγ·f(x)in the weight-sum approach can be interpreted as the attempt toﬁnd the y value for which,starting from the origin point,the line with slope−γ1γ2is tangent to the region C.Obviously,changing the weight parameters leads to possibly different touching points of the line to the feasible region.If the Pareto curve is convex then there is room to calculate such points for differentγvectors(see Fig.2.3).2 f1(xFig.2.3Geometrical representation of the weight-sum approach in the convex Pareto curve caseOn the contrary,when the curve is non-convex,there is a set of points that cannot be reached for any combinations of theγweight vector(see Fig.2.4).182Multi-objective Optimization f1(xFig.2.4Geometrical representation of the weight-sum approach in the non-convex Pareto curve caseThe following result by Geoffrion(1968)states a necessary and sufﬁcient condi-tion in the case of convexity as follows:If the solution set S is convex and the n objectives f i are convex on S,x∗is a strict Pareto optimum if and only if it existsγ∈R n,such that x∗is an optimal solution of problem P s(γ).Similarly:If the solution set S is convex and the n objectives f i are convex on S,x∗is a weak Pareto optimum if and only if it existsγ∈R n,such that x∗is an optimal solution of problem P w(γ).If the convexity hypothesis does not hold,then only the necessary condition re-mains valid,i.e.,the optimal solutions of P s(γ)and P w(γ)are strict and weak Pareto optima,respectively.2.3.2ε-constraints MethodBesides the scalarization approach,another solution technique to multi-objective optimization is theε-constraints method proposed by Chankong and Haimes in 1983.Here,the decision maker chooses one objective out of n to be minimized; the remaining objectives are constrained to be less than or equal to given target val-2.3Techniques to Solve Multi-objective Optimization Problems19 ues.In mathematical terms,if we let f2(x)be the objective function chosen to be minimized,we have the following problem P(ε2):min f2(x)f i(x)≤εi,∀i∈{1,...,n}\{2}x∈S.We note that this formulation of theε-constraints method can be derived by a more general result by Miettinen,that in1994proved that:If an objective j and a vectorε=(ε1,...,εj−1,εj+1,...,εn)∈R n−1exist,such that x∗is an optimal solution to the following problem P(ε):min f j(x)f i(x)≤εi,∀i∈{1,...,n}\{j}x∈S,then x∗is a weak Pareto optimum.In turn,the Miettinen theorem derives from a more general theorem by Yu(1974) stating that:x∗is a strict Pareto optimum if and only if for each objective j,with j=1,...,n, there exists a vectorε=(ε1,...,εj−1,εj+1,...,εn)∈R n−1such that f(x∗)is the unique objective vector corresponding to the optimal solution to problem P(ε).Note that the Miettinen theorem is an easy implementable version of the result by Yu(1974).Indeed,one of the difﬁculties of the result by Yu,stems from the uniqueness constraint.The weaker result by Miettinen allows one to use a necessary condition to calculate weak Pareto optima independently from the uniqueness of the optimal solutions.However,if the set S and the objectives are convex this result becomes a necessary and sufﬁcient condition for weak Pareto optima.When,as in problem P(ε2),the objective isﬁxed,on the one hand,we have a more simpliﬁed version,and therefore a version that can be more easily implemented in automated decision-support systems;on the other hand,however,we cannot say that in the presence of S convex and f i convex,∀i=1,...,n,all the set of weak Pareto optima can be calculated by varying theεvector.One advantage of theε-constraints method is that it is able to achieve efﬁcient points in a non-convex Pareto curve.For instance,assume we have two objective202Multi-objective Optimization functions where objective function f1(x)is chosen to be minimized,i.e.,the problem ismin f1(x)f2(x)≤ε2x∈S,we can be in the situation depicted in Fig.2.5where,when f2(x)=ε2,f1(x)is an efﬁcient point of the non-convex Pareto curve.f1(xf 2(x)£e2x)f1(xFig.2.5Geometrical representation of theε-constraints approach in the non-convex Pareto curve caseTherefore,as proposed in Steurer(1986)the decision maker can vary the upper boundsεi to obtain weak Pareto optima.Clearly,this is also a drawback of this method,i.e.,the decision maker has to choose appropriate upper bounds for the constraints,i.e.,theεi values.Moreover,the method is not particularly efﬁcient if the number of the objective functions is greater than two.For these reasons,Erghott and Rusika in2005,proposed two modiﬁcations to improve this method,with particular attention to the computational difﬁculties that the method generates.2.3Techniques to Solve Multi-objective Optimization Problems21 2.3.3Goal ProgrammingGoal Programming dates back to Charnes et al.(1955)and Charnes and Cooper (1961).It does not pose the question of maximizing multiple objectives,but rather it attempts toﬁnd speciﬁc goal values of these objectives.An example can be given by the following program:f1(x)≥v1f2(x)=v2f3(x)≤v3x∈S.Clearly we have to distinguish two cases,i.e.,if the intersection between the image set C and the utopian set,i.e.,the image of the admissible solutions for the objectives,is empty or not.In the former case,the problem transforms into one in which we have toﬁnd a solution whose value is as close as possible to the utopian set.To do this,additional variables and constraints are introduced.In particular,for each constraint of the typef1(x)≥v1we introduce a variable s−1such that the above constraint becomesf1(x)+s−1≥v1.For each constraint of the typef2(x)=v2we introduce a surplus two variables s+2and s−2such that the above constraint be-comesf2(x)+s−2−s+2=v2.For each constraint of the typef3(x)≤v3we introduce a variable s+3such that the above constraint becomesf3(x)−s+3≤v3.222Multi-objective OptimizationLet us denote with s the vector of the additional variables.A solution(x,s)to the above problem is called a strict Pareto-slack optimum if and only if a solution (x ,s ),for every x ∈S,such that s i≤s i with at least one strict inequality does not exist.There are different ways of optimizing the slack/surplus variables.An exam-ple is given by the Archimedean goal programming,where the problem becomes that of minimizing a linear combination of the surplus and slack variables each one weighted by a positive coefﬁcientαas follows:minαs−1s−1+αs+2s+2+αs−2s−2+αs+3s+3f1(x)+s−1≥v1f2(x)+s−2−s+2=v2f3(x)−s+3≤v3s−1≥0s+2≥0s−2≥0s+3≥0x∈S.For the above problem,the Geoffrion theorem says that the resolution of this prob-lem offers strict or weak Pareto-slack optimum.Besides Archimedean goal programming,other approaches are the lexicograph-ical goal programming,the interactive goal programming,the reference goal pro-gramming and the multi-criteria goal programming(see,e.g.,T’kindt and Billaut, 2005).2.3.4Multi-level ProgrammingMulti-level programming is another approach to multi-objective optimization and aims toﬁnd one optimal point in the entire Pareto surface.Multi-level programming orders the n objectives according to a hierarchy.Firstly,the minimizers of theﬁrst objective function are found;secondly,the minimizers of the second most important2.3Techniques to Solve Multi-objective Optimization Problems23objective are searched for,and so forth until all the objective function have been optimized on successively smaller sets.Multi-level programming is a useful approach if the hierarchical order among the objectives is meaningful and the user is not interested in the continuous trade-off among the functions.One drawback is that optimization problems that are solved near the end of the hierarchy can be largely constrained and could become infeasi-ble,meaning that the less important objective functions tend to have no inﬂuence on the overall optimal solution.Bi-level programming(see,e.g.,Bialas and Karwan,1984)is the scenario in which n=2and has received several attention,also for the numerous applications in which it is involved.An example is given by hazmat transportation in which it has been mainly used to model the network design problem considering the government and the carriers points of view:see,e.g.,the papers of Kara and Verter(2004),and of Erkut and Gzara(2008)for two applications(see also Chap.4of this book).In a bi-level mathematical program one is concerned with two optimization prob-lems where the feasible region of theﬁrst problem,called the upper-level(or leader) problem,is determined by the knowledge of the other optimization problem,called the lower-level(or follower)problem.Problems that naturally can be modelled by means of bi-level programming are those for which variables of theﬁrst problem are constrained to be the optimal solution of the lower-level problem.In general,bi-level optimization is issued to cope with problems with two deci-sion makers in which the optimal decision of one of them(the leader)is constrained by the decision of the second decision maker(the follower).The second-level de-cision maker optimizes his/her objective function under a feasible region that is deﬁned by theﬁrst-level decision maker.The latter,with this setting,is in charge to deﬁne all the possible reactions of the second-level decision maker and selects those values for the variable controlled by the follower that produce the best outcome for his/her objective function.A bi-level program can be formulated as follows:min f(x1,x2)x1∈X1x2∈argmin g(x1,x2)x2∈X2.242Multi-objective OptimizationThe analyst should pay particular attention when using bi-level optimization(or multi-level optimization in general)in studying the uniqueness of the solutions of the follower problem.Assume,for instance,one has to calculate an optimal solu-tion x∗1to the leader model.Let x∗2be an optimal solution of the follower problem associated with x∗1.If x∗2is not unique,i.e.,|argmin g(x∗1,x2)|>1,we can have a sit-uation in which the follower decision maker can be free,without violating the leader constraints,to adopt for his problem another optimal solution different from x∗2,i.e.,ˆx2∈argmin g(x∗1,x2)withˆx2=x∗2,possibly inducing a f(x∗1,ˆx2)>f(x∗1,x∗2)on the leader,forcing the latter to carry out a sensitivity analysis on the values at-tained by his objective function in correspondence to all the optimal solutions in argmin g(x∗1,x2).Bi-level programs are very closely related to the van Stackelberg equilibrium problem(van Stackelberg,1952),and the mathematical programs with equilibrium constraints(see,e.g.,Luo et al.1996).The most studied instances of bi-level pro-gramming problems have been for a long time the linear bi-level programs,and therefore this subclass is the subject of several dedicated surveys,such as that by Wen and Hsu(1991).Over the years,more complex bi-level programs were studied and even those including discrete variables received some attention,see,e.g.,Vicente et al.(1996). Hence,more general surveys appeared,such as those by Vicente and Calamai(1994) and Falk and Liu(1995)on non-linear bi-level programming.The combinatorial nature of bi-level programming has been reviewed in Marcotte and Savard(2005).Bi-level programs are hard to solve.In particular,linear bi-level programming has been proved to be strongly NP-hard(see,Hansen et al.,1992);Vicente et al. (1996)strengthened this result by showing thatﬁnding a certiﬁcate of local opti-mality is also strongly NP-hard.Existing methods for bi-level programs can be distinguished into two classes.On the one hand,we have convergent algorithms for general bi-level programs with the-oretical properties guaranteeing suitable stationary conditions;see,e.g.,the implicit function approach by Outrata et al.(1998),the quadratic one-level reformulation by Scholtes and Stohr(1999),and the smoothing approaches by Fukushima and Pang (1999)and Dussault et al.(2004).With respect to the optimization problems with complementarity constraints, which represent a special way of solving bi-level programs,we can mention the pa-pers of Kocvara and Outrata(2004),Bouza and Still(2007),and Lin and Fukushima2.4Multi-objective Optimization Integer Problems25(2003,2005).Theﬁrst work presents a new theoretical framework with the im-plicit programming approach.The second one studies convergence properties of a smoothing method that allows the characterization of local minimizers where all the functions deﬁning the model are twice differentiable.Finally,Lin and Fukushima (2003,2005)present two relaxation methods.Exact algorithms have been proposed for special classes of bi-level programs, e.g.,see the vertex enumeration methods by Candler and Townsley(1982),Bialas and Karwan(1984),and Tuy et al.(1993)applied when the property of an extremal solution in bi-level linear program plementary pivoting approaches(see, e.g.,Bialas et al.,1980,and J´u dice and Faustino,1992)have been proposed on the single-level optimization problem obtained by replacing the second-level optimiza-tion problem by its optimality conditions.Exploiting the complementarity structure of this single-level reformulation,Bard and Moore(1990)and Hansen et al.(1992), have proposed branch-and-bound algorithms that appear to be among the most efﬁ-cient.Typically,branch-and-bound is used when the lower-level problem is convex and regular,since the latter can be replaced by its Karush–Kuhn–Tucker(KKT) conditions,yielding a single-level reformulation.When one deals with linear bi-level programs,the complementarity conditions are intrinsically combinatorial,and in such cases branch-and-bound is the best approach to solve this problem(see,e.g., Colson et al.,2005).A cutting-plane approach is not frequently used to solve bi-level linear programs.Cutting-plane methods found in the literature are essentially based on Tuy’s concavity cuts(Tuy,1964).White and Anandalingam(1993)use these cuts in a penalty function approach for solving bi-level linear programs.Marcotte et al.(1993)propose a cutting-plane algorithm for solving bi-level linear programs with a guarantee ofﬁnite termination.Recently,Audet et al.(2007),exploiting the equivalence of the latter problem with a mixed integer linear programming one, proposed a new branch-and-bound algorithm embedding Gomory cuts for bi-level linear programming.2.4Multi-objective Optimization Integer ProblemsIn the previous section,we gave general results for continuous multi-objective prob-lems.In this section,we focus our attention on what happens if the optimization problem being solved has integrality constraints on the variables.In particular,all262Multi-objective Optimizationthe techniques presented can be applied in these situations as well,with some lim-itations on the capabilities of these methods to construct the Pareto front entirely. Indeed,these methods are,in general,very hard to solve in real applications,or are unable toﬁnd all efﬁcient solutions.When integrality constraints arise,one of the main limits of these techniques is in the inability of obtaining some Pareto optima; therefore,we will have supported and unsupported Pareto optima.f 2(x)f1(xFig.2.6Supported and unsupported Pareto optimaFig.2.6gives an example of these situations:points p6and p7are unsupported Pareto optima,while p1and p5are supported weak Pareto optima,and p2,p3,and p4are supported strict Pareto optima.Given a multi-objective optimization integer problem(MOIP),the scalarization in a single objective problem with additional variables and/or parameters toﬁnd a subset of efﬁcient solutions to the original MOIP,has the same computational complexity issues of a continuous scalarized problem.In the2006paper of Ehrgott“A discussion of scalarization techniques for mul-tiple objective integer programming”the author,besides the scalarization tech-niques also presented in the previous section(e.g.,the weighted-sum method,the ε-constraint method),satisfying the linear requirement imposed by the MOIP for-mulation(where variables are integers,but constraints and objectives are linear),2.4Multi-objective Optimization Integer Problems27presented more methods like the Lagrangian relaxation and the elastic-constraints method.By the author’s analysis,it emerges that the attempt to solve the scalarized prob-lem by means of Lagrangian relaxation would not lead to results that go beyond the performance of the weighted-sum technique.It is also shown that the general linear scalarization formulation is NP-hard.Then,the author presents the elastic-constraints method,a new scalarization technique able to overcome the drawback of the previously mentioned techniques related toﬁnding all efﬁcient solutions,com-bining the advantages of the weighted-sum and theε-constraint methods.Further-more,it is shown that a proper application of this method can also give reasonable computing times in practical applications;indeed,the results obtained by the author on the elastic-constraints method are applied to an airline-crew scheduling problem, whose size oscillates from500to2000constraints,showing the effectiveness of the proposed technique.2.4.1Multi-objective Shortest PathsGiven a directed graph G=(V,A),an origin s∈V and a destination t∈V,the shortest-path problem(SPP)aims toﬁnd the minimum distance path in G from o to d.This problem has been studied for more than50years,and several polynomial algorithms have been produced(see,for instance,Cormen et al.,2001).From the freight distribution point of view the term shortest may have quite dif-ferent meanings from faster,to quickest,to safest,and so on,focusing the attention on what the labels of the arc set A represent to the decision maker.For this reason, in some cases we willﬁnd it simpler to deﬁne for each arc more labels so as to represent the different arc features(e.g.,length,travel time,estimated risk).The problem toﬁnd multi-objective shortest paths(MOSPP)is known to be NP-hard(see,e.g.,Seraﬁni,1986),and the algorithms proposed in the literature faced the difﬁculty to manage the large number of non-dominated paths that results in a considerable computational time,even in the case of small instances.Note that the number of non-dominated paths may increase exponentially with the number of nodes in the graph(Hansen,1979).In the multi-objective scenario,each arc(i,j)in the graph has a vector of costs c i j∈R n with c i j=(c1i j,...,c n i j)components,where n is the number of criteria.。

Abstract

பைடு நூலகம்
447
argument has been understood to undergo the syntactic/pragmatic ellipsis. However, we observe some delicate differences in intuition on whether or not the underlined constituents in (1-6) can be omitted:
(2)
(3)t
(4)
(5)
(6)
Try to understand the above sentences, omitting the underlined constituents. If someone omits °1 i dohoeji-e 'in this urban area' in sentence (1), the resulting sentence probably makes us to be anxious to know where the event happened, though we do not think that the sentence is wrong (i.e. ungrammatical). Example (with the underlined part omitted) (2) might lead us to ask "When? or Why?", though it is grammatical. As for example (3), we are likely to use it pragmatically, that is, supposing that we already know where he stayed. Example (4a) is similar to (3). It is more natural if it is colloquially uttered like (4b). Example (5b) seems to be more natural than (5a), though both have the same verb 11114 byenhaessda 'turned'. Readers are likely to understand (5b), guessing diverse situations such as Cheolsu got a cold, he was at the age of the voice change, or he got hoarse because he used vocal chords too much. In contrast, example (5a) is very unnatural without an appropriate context. Example (6) seems to be natural because we can obviously infer the omitted constituent saegeos-eulo 'into a new one'. One may not even recognize the ellipsis itself in (6). In distinguishing optional arguments from obligatory ones, Ki-Shim Nam (1993) notes the following: /o 'to' is obligatory or not in a given "In many cases, it is not easy to decide whether sentence. There are no formal criteria. There is only a heuristic method such as: The given /o 'to' is obligatory in a sentence if and only if the elliptical structure becomes ungrammatical. To decide the possibility of its ellipsis in a sentence, we do nothing but consider the semantic property of the predicate, totally depending on our linguistic intuition." The problem is that native speakers are split into various groups on whether or not the examples in (1-6), without the underlined parts are right (i.e. natural or grammatical): total agreement, a little bit

拓扑优化在结构工程中的应用

第36卷第6期2020年12月结构工程师Structural Engineers Vol.36，No.6Dec.2020拓扑优化在结构工程中的应用高文俊1，2，*吕西林1，2（1.同济大学土木工程防灾国家重点实验室，上海200092；2.同济大学结构防灾减灾工程系，上海200092）摘要拓扑优化技术经过多年的发展已成为结构设计的有力工具。

在过去的十年中，拓扑优化在结构工程领域内涌现出一批具有创新性的应用。

从结构理论到构件设计，再到整体结构找形，这些应用涉及工程结构的各个层面。

拓扑优化在这些应用中被视为一种突破传统设计的重要方法。

本文对拓扑优化在结构工程中的应用进行了归纳与梳理，并对拓扑优化方法在钢筋混凝土构件设计、结构构件布置、构件形态设计、结构找形中具有代表性的应用研究进行了介绍；然后对其发展现状进行了讨论；最后对其发展趋势与潜在的应用方向进行了展望。

关键词拓扑优化，结构找形，结构工程，工程应用Applications of Topology Optimization in Structural EngineeringGAO Wenjun1，2LÜXilin1，2（1.State Key Laboratory of Disaster Reduction in Civil Engineering，Shanghai200092，China；2.Department of DisasterMitigation for Structures，Tongji University，Shanghai200092，China）Abstract Topology optimization has undergone a long period of development and become a powerful tool for structural design.Over the last decade，topology optimization inspired innovative applications in structural engineering.These applications involve multiple levels of engineering structures，from structural theories to designs of structural members，then to form finding of a whole structure.Topology optimization is viewed as an important method for breaking through traditional designs in these applications.This paper summarizes these applications，and also provides an overview of some representative ones in designing reinforced concrete members，layout of structural members，morphology of structural components and structural form finding.The current status of applying topology optimization in structural engineering is discussed.The development tendency and potential applications are prospected at the end of this paper.Keywords topology optimization，form finding of structures，structural engineering，applications in engineering0引言为了实现结构创新以满足新的社会需求，从20世纪70年代开始，结构工程不断吸收与融合工程力学、计算机科学、材料学、机电控制、信息工程、智能技术等学科的先进技术，从而极大地丰富了自身的内涵并取得了许多创新成果。

SIAMJ.Optimization DISCRETIZEDCONTINUOUSMINIMAXPROBLEMSANDOTHER Vol.x,No.x,pp.x{xx,xxx19xx

H = fH = H T : 1 kdk2 hd; Hdi
2 kdk2
8d 2 IRn g:
let d(x; H; ^ ) be the unique KKT point of QP (x; H; ^ ). Then Lemma 3.2 implies that d(x ; H; ^ ) = 0 for all H 2 H. Since H is compact, in view of Assumptions 4 and 5, it follows from a classical result of Robinson's 1, Theorem 2.1] that, given > 0, there exists ^ > 0 such that, for all x with kx ? x k < ^ and all H 2 H, all ! 2 max (x ) are binding for QP (x; H; ^ ) and kd(x; H; ^ )k < . That ^ > 0 can be chosen independent of ^ follows from niteness of . b. Lemma 3.14'. For k large enough, max(x ) k Proof. Let 1 ; 2 > 0 be as given by Assumption 3, and let > 0 be as given by Lemma 3.13' (for an arbitrary > 0). Since xk ! x as k ! 1 (Proposition 3.11) ^ ^ it follows from Lemma 3.12 that there exists k such that kxk ? x k < for all k k

梯度下降算法研究综述

文章编号：2096-1472(2020)-02-01-04DOI:10.19644/ki.issn2096-1472.2020.02.001软件工程 SOFTWARE ENGINEERING 第23卷第2期2020年2月V ol.23 No.2Feb. 2020梯度下降算法研究综述李兴怡，岳洋(上海理工大学，上海 200093)摘要：在机器学习领域中，梯度下降算法是一种广泛用于求解线性和非线性模型最优解的迭代算法，它的中心思想在于通过迭代次数的递增，调整使得损失函数最小化的权重。

本文首先概述了基于多元线性模型的梯度下降算法；其次介绍了梯度下降算法三种框架，使用Python实现了自主停止训练的BGD算法；针对梯度下降算法存在的不足，综述了近三年算法优化的研究成果。

最后，总结了本文的主要研究工作,对梯度下降优化算法的研究趋势进行了展望。

关键词：机器学习；多元线性模型；梯度下降算法；算法实现；优化算法中图分类号：TP181 文献标识码：ASurvey of Gradient Descent AlgorithmLI Xingyi,YUE Yang(University of Shanghai for Science and Technology ,Shanghai 200093,China )Abstract:Gradient descent algorithm is an iteration algorithm which is widely used in figuring out the minimum of linear model and non-linear model in the field of machine learning.Its main idea is to adjust the weights that minimizes the cost function through increasing the number of iterations.First,this paper outlines gradient descent algorithm based on multivariate linear model.Then introduces three kinds of gradient descent variants and accomplishes batch gradient descent algorithms which stops training autonomously by means of Python.This paper also gives an overview of research advances from recent three years of gradient descent optimization algorithms to existing deficiencies.At the end of this paper,the main research findings are summarized and the research tendency of gradient descent optimization algorithms is prospected.Keywords:machine learning;multivariate linear model;gradient descent algorithm;algorithm implementation; optimization algorithm1 引言(Introduction)在求解机器学习中无约束优化问题的方法中，优化方法是必不可少的一环。

萝卜籽酶解制备萝卜硫素工艺优化及其体外消化研究

李惠琳，刘昊，李珏，等. 萝卜籽酶解制备萝卜硫素工艺优化及其体外消化研究[J]. 食品工业科技，2024，45（9）：159−167. doi:10.13386/j.issn1002-0306.2023050037LI Huilin, LIU Hao, LI Jue, et al. Process Optimization and in Vitro Digestion Research of Raphanus sativus Seeds Sulforaphane Prepared by Enzymolysis Method[J]. Science and Technology of Food Industry, 2024, 45(9): 159−167. (in Chinese with English abstract). doi: 10.13386/j.issn1002-0306.2023050037· 工艺技术 ·萝卜籽酶解制备萝卜硫素工艺优化及其体外消化研究李惠琳1,2，刘　昊1，李　珏1，木耶赛尔•凯代斯1，图尔荪阿依•图尔贡1，赵　雷1，何庆峰1,2,*（1.和田职业技术学院农业科技系，新疆维吾尔自治区和田 848000；2.天津农学院食品科学与生物工程学院，天津 300000）摘　要：为探究酶解法制备萝卜籽萝卜硫素的最佳工艺及其消化特性。

采用红萝卜种子为原料，以萝卜硫素得率为指标，通过单因素实验和响应面优化试验得出最佳酶解工艺条件。

进一步采用体外模拟胃肠消化模型，探究萝卜籽乙酸乙酯提取物萝卜硫素含量及其抗氧化活性的变化规律。

得到萝卜籽制备萝卜硫素的最佳酶解条件为：酶解时间22 min ，酶解温度40 ℃，V C 添加量0.8 mg/g ，在此最佳工艺条件下，萝卜硫素得率为2.11±0.02 mg/g ，与预测值（2.14 mg/g ）接近，相对误差为1.4%，酶解工艺切实可行。

响应面法优化准噶尔山楂总三萜提取工艺及其纯化工艺研究

冯琳，罗世博，古丽格娜·皮达买买提，等. 响应面法优化准噶尔山楂总三萜提取工艺及其纯化工艺研究[J]. 食品工业科技，2024，45（9）：196−204. doi: 10.13386/j.issn1002-0306.2023070139FENG Lin, LUO Shibo, PIDAMAIMAITI Guligena, et al. Optimization of Extraction Process of Total Triterpenoids from Crataegus songarica by Response Surface Methodology and Its Purification Process[J]. Science and Technology of Food Industry, 2024, 45(9):196−204. (in Chinese with English abstract). doi: 10.13386/j.issn1002-0306.2023070139· 工艺技术 ·响应面法优化准噶尔山楂总三萜提取工艺及其纯化工艺研究冯　琳1,2，罗世博1,2，古丽格娜·皮达买买提1，刘媛梦1,3，高红艳1,2,*（1.伊犁师范大学化学化工学院，新疆伊宁 835000；2.新疆普通高等学校天然产物化学与应用重点实验室，新疆伊宁 835000；3.污染物化学与环境治理重点实验室，新疆伊宁 835000）摘　要：为了提高准噶尔山楂的高值化利用，以准噶尔山楂为原料，通过单因素实验、正交试验和响应面设计优化复合酶辅助提取准噶尔山楂总三萜的工艺，并对纯化工艺进行优化。

结果表明，准噶尔山楂总三萜的最优提取工艺为纤维素酶、果胶酶和木瓜蛋白酶的添加量分别为4%、4%和2%，料液比为1:24 g/mL ，温度51 ℃，pH 为4.5，酶解76 min ，此条件下准噶尔山楂总三萜的提取量为（36.570±0.332）mg/g ，通过超声再次提取，总三萜的提取量达到（53.782±0.673）mg/g 。

Optimization Toolbox MATLAB优化工具箱

Optimization Toolbox－－求解常规和大型优化问题Optimization Toolbox 提供了应用广泛的算法集合，用于求解常规和大型的优化问题。

这些算法解决带约束、无约束的、连续的和离散的优化问题。

这些算法可以求解带约束的、无约束的以及离散的优化问题。

工具箱中包含的函数可以用于线性规划、二次规划、二进制整数规划、非线性优化、非线性最小二乘、非线性方程、以及多目标优化等。

用户能够使用这些算法寻找最优解，进行权衡分析，在多个设计方案之间平衡，以及将优化算法集成到算法和模型之中。

主要特点∙交互式工具用于定义、求解优化问题，并能监控求解过程∙求解非线性优化和多目标优化问题∙求解非线性最小二乘，数据拟合和非线性方程∙提供了解决二次方程和线性规划问题的方法∙提供了解决二进制整数规划问题的方法∙某些带约束条件的非线性求解器支持并行运算使用Optimization Toolbox 中的基于梯度的求解器寻找峰值函数（peaks function）的局部最小解。

运用优化工具箱提供的大型线性最小二乘法修复一张模糊的照片。

定义，求解以及评定优化问题优化工具箱提供了解决极小极大值问题的最常用方法。

工具箱包含了常规和大型优化问题的算法，使用户可以利用问题的稀疏结构来求解问题。

用户可以通过命令行或图形用户界面Optimization Tool调用工具箱函数和求解器选项。

通过命令行运行的优化程序(左，调用了定义指标函数（右上）和限定条件方程（右下）的MATLAB文件。

Optimization Tool 是一个将一般优化工作简单化的图形用户界面。

通过该图形用户界面，用户能够完成以下操作：∙定义自己的优化问题并选择求解器∙配置，检验优化选项和所选求解器的默认设置∙运行优化问题，显示中间以及最终结果∙在可选择的快速帮助窗口中查看特定求解器的文档∙在MATLAB 的工作空间和优化工具之间导入和导出用户问题的定义，算法配置和结果∙保存用户工作和使工作自动化，自动生成M 语言代码∙调用Global Optimization Toolbox中的求解器使用Optimization Tool 设置并求解的一个优化程序(左)。

1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。

Discrete OptimizationAn alternative framework to Lagrangian relaxationapproach for job shop schedulingHaoxun Chen *,Peter B.LuhDepartment of Electrical and Computer Engineering,University of Connecticut,Storrs,CT 06269,USAReceived 9May 2000;accepted 8May 2002AbstractA new Lagrangian relaxation (LR)approach is developed for job shop scheduling problems.In the approach,operation precedence constraints rather than machine capacity constraints are relaxed.The relaxed problem is de-composed into single or parallel machine scheduling subproblems.These subproblems,which are NP-complete in general,are approximately solved by using fast heuristic algorithms.The dual problem is solved by using a recently developed ‘‘surrogate subgradient method’’that allows approximate optimization of the subproblems.Since the al-gorithms for subproblems do not depend on the time horizon of the scheduling problems and are very fast,our new LR approach is eﬃcient,particularly for large problems with long time horizons.For these problems,the machine de-composition-based LR approach requires much less memory and computation time as compared to a part decom-position-based approach as demonstrated by numerical testing.Ó2002Elsevier Science B.V.All rights reserved.Keywords:Scheduling;Job shops;Lagrangian relaxation;Manufacturing systems1.IntroductionScheduling is to allocate ﬁnite resources over time to accomplish a given set of tasks.High quality schedules can improve delivery perfor-mance,reduce inventory costs,and are very im-portant to manufacturers in today Õs time-based competition.To obtain a high quality schedulewithin an acceptable computation time,however,is extremely diﬃcult because of the NP-hard na-ture of the problem and the large sizes of practical applications [12].Scheduling methods in the literature can be classiﬁed into optimization methods,approximate and heuristics methods,and dispatching rules.Optimization methods include branch and bound and dynamic programming (DP)(e.g.,[3]and [17]).Although they can obtain optimal schedules,they are computationally intensive in general.By contrast,approximate and heuristics methods make a tradeoﬀbetween solution quality and computation time.These methods include simu-lated annealing [11],genetic algorithms [7,8,13],*Corresponding author.Current Address:Laboratoired ÕOptimisation des Syste`mes Industriels,Universite ´de Tech-nologie de Troyes,12Rue Marie Curie,BP 20,10010Troyes Cedex,France.Tel.:+33-325-715642;fax:+33-325-715649.E-mail addresses:hxchen@,chen@utt.fr (H.Chen),luh@ (P.B.Luh).0377-2217/02/$-see front matter Ó2002Elsevier Science B.V.All rights reserved.doi:10.1016/S0377-2217(02)00470-8European Journal of Operational Research 149(2003)499–512/locate/dswtabu search[6,14],constraint programming[15], shifting bottleneck procedure[1],and Lagrangian relaxation[4,10,19].Dispatching rules are too simple to be eﬀective in most cases[2,16].Lagrangian relaxation(LR)has recently emer-ged as a practical approach for complex schedul-ing problems[4,10,19].It is a decomposition and coordination approach which can obtain near optimal schedules with quantiﬁable quality in a reasonable computation time for practical sched-uling problems.n this approach,machine capac-ity constraints areﬁrst relaxed by using Lagrange multipliers.The relaxed problem can be decom-posed into a set of part subproblems that are solved by using DP.The multipliers are then iter-atively adjusted at the‘‘high level’’based on the degree of constraint violation.At the termination of such iterations,a simple heuristic is applied to adjust subproblem solutions to obtain a feasible schedule satisfying all constraints.This approach is quite eﬃcient in general. However,as the computation complexity of the DP algorithm and the number of multipliers are proportional to the time horizon,the approach becomes time and memory consuming for prob-lems with long time horizons.This paper presents a new LR approach for job shop scheduling problems to minimize a weighted earliness and tardiness criterion.Our goal is toﬁnd good schedules for problems of large size in a short computation time.First,a separable formulation of the problems is presented in Section2.The new approach is then introduced in Section3.In the approach,operation precedence constraints rather than machine capacity constraints are relaxed.The relaxed problem is decomposed into single or parallel machine scheduling subproblems.These subproblems are NP-complete in general,and are approximately solved by using fast heuristic algo-rithms as presented in Section4.The dual problem is iteratively solved by using a recently developed ‘‘surrogate subgradient(SSG)method’’that allows approximate optimization of the subprob-lems[20].Feasible schedules are then constructed by using heuristics based on subproblem solutions or multipliers.The solution of the dual problem and the construction of feasible schedules will be presented in Section5.Because the subproblem algorithms do not de-pend on the time horizon and are very fast,our new LR approach is eﬃcient,particularly for problems with long time horizons.For these problems,the new approach requires much less memory and computation time as compared to the LR approach based on part decomposition.This is demonstrated by numerical testing presented in Section6.Test-ing also shows that our new approach signiﬁcantly outperforms dispatching rules and can solve problems with thousands of parts and tens of ma-chine types in less than15minutes on a personal computer,making it practical for shop-ﬂoor use.2.Problem formulationThe job shop scheduling problem considered in this paper is to schedule N parts on H types of machines to minimize a weighted earliness and tardiness criterion following[19].Each machine type hð16h6HÞhas M h identical machines,and the completion of each part ið16i6NÞrequires a series of N i operations,denoted byði;1Þ;ði;2Þ;...;ði;N iÞ.Each operationði;jÞis non-preemptive and can be performed on a machine belonging to a set of alternative machine types H ij f1;2;...;H g. For simplicity,all parts are assumed to be avail-able at time0.For the case where some parts are not available at time0,the following formulation and solution methodology are still applicable after a slight modiﬁcation.The following symbols will be used in the problem formulation:B i beginning time of part ic ij completion time of operationði;jÞC i completion time of part iD i due date of part iE i earliness of part i,deﬁned as maxð0;D iÀC iÞm ij machine type selected to process opera-tionði;jÞ,m ij2H ijM h s number of available machines of type h at time s,16h6HO h set of operations that can be performed on machine type ht ijh processing time of operationði;jÞon ma-chine type h2H ij500H.Chen,P.B.Luh/European Journal of Operational Research149(2003)499–512T i tardiness of part i,deﬁned as maxð0;C iÀD iÞw i tardiness weight for part ib i earliness weight for part id ijh operation index over time for machinetype h with d ijhðsÞ¼1if operationði;jÞisperformed on a machine of type h attime s,and d ijhðsÞ¼0otherwise.That is,d ijhðsÞ¼1forði;j;h;sÞwith m ij¼h andc ijÀt ijh6s6c ij,and d ijhðsÞ¼0for othersWith the above symbols,an optimization model for the scheduling problem is presented below following[19]:Objective.The goal of on-time delivery and low work-in-process inventory is modeled as a weigh-ted tardiness and earliness cost,i.e.,min f m ij g;f c ij g J;with J¼Xiw i T iþXibiE i;ð1Þwhere w i and b i are tardiness and earliness weights, T i and E i are tardiness and earliness for part i, respectively.Since on-time delivery is the foremost criterion in(1),b i is usually an order of magnitude smaller than w i.Operation precedence constraints.An operation cannot be started until its preceding operation is ﬁnished,and it requires a speciﬁc amount of time for processing on the selected machine type,i.e.,c i;jÀ1þt ijmij 6ci;j;i¼1;2;...;N;j¼1;2;...;N i;ð2Þwhere c i;j is the completion time of operationði;jÞ,t ijmij is the processing time of operationði;jÞonmachine type m ij.Machine capacity constraints.The number of operations being processed on a machine of type h at any time instant s cannot exceed M h,the number of machines for the type,i.e.,Xijd ijhðsÞ6M h s;06s<1;h¼1;2;...;H;ð3Þwhere d ijh is operation index over time for machine type h,M h s is the number of available machines of type h at time s.The above formulation diﬀers from that of[19] in its time concept,as continuous time rather than discrete time is used.Consequently,the left side of the inequality(2)does not contain a term)1.The overall problem is to minimize the cost function(1)subject to the above constraints by selecting appropriate machine types and comple-tion times for individual operations.Since the problem is NP-hard and no algorithm can opti-mally solve the problem of practical sizes in a reasonable computation time,a near-optimal ap-proach based on Lagrangian relaxation will be developed in the following sections.3.Solution frameworkOur new LR approach is based on machine decomposition,and is carried out by relaxing op-eration precedence constraints.The reason for relaxing the precedence constraints is that the number of such constraints does not depend on the time horizon so that relaxation can lead to a time horizon-independent approach.However,if only precedence constraints are relaxed,the re-laxed subproblems,which are parallel machine problems to minimize the weighted earliness and tardiness criterion,are diﬃcult to solve.For this reason,an additional eﬀort is made by reformu-lating the original problem and relaxing also the earliness and tardiness constraints.The resulting subproblems are then parallel machine problems to minimize the weighted completion time,where eﬃcient approximate algorithms exist.3.1.Problem reformulationThe deﬁnition of the tardiness for part i implies thatC iÀD i6T i;andð4aÞ06T i;i¼1;2;...;N:ð4bÞSimilarly,the deﬁnition of the earliness for part i implies thatD iÀC i6E i;ð5aÞ06E i;i¼1;2;...;N:ð5bÞH.Chen,P.B.Luh/European Journal of Operational Research149(2003)499–512501With the above derivation,the problem presented in the previous section can be reformulated as follows:P:minf m ij g;f c ij g;f T i g;f E ig J;with J¼Xiw i T iþXib i E i;subject to(2),(3),(4a),(4b),(5a)and(5b),where earliness E i and tardiness T i are also treated as additional decision variables.These two problems are equivalent because for any optimal solution of P,T i¼maxð0;C iÀD iÞand E i¼maxð0;D iÀC iÞfor all i.If this is not true, there exists some i such that T i>maxð0;C iÀD iÞor E i>maxð0;D iÀC iÞ.A better solution can then be obtained by taking T i¼maxð0;C iÀD iÞandE i¼maxð0;D iÀC iÞ.3.2.Relaxation and decompositionBy introducing multipliers f l i g,f m i g,and f k ij g to relax constraints(4a),(5a),and(2),respectively, the relaxed problem of P is obtained as:RP:minf m ij g;f c ij g;f T i g;f E i gLðff l i g;f m i g;f k ij gg;ff m ij g;f c ij g;f T i g;f E i ggÞ;ð6Þsubject to(4b),(5b),(3),withL JþXi liðC iÀD iÀT iÞþXim iðD iÀC iÀE iÞþXij k ijðc i;jÀ1þt ijmijÀc i;jÞ:Let ff mÃij g;f cÃijg;f TÃig;f EÃigg be an optimal solu-tion of RP for a given set of multipliers.The dual problem of RP is to maximize the dual function q: DP:maxf l i;m i;k ijg P0qð7Þwith q Lðff l i g;f m i g;f k ij gg;ff mÃij g;f cÃijg;f TÃig;f EÃi ggÞ.Function L can be rewritten as:L¼Xi ðw iÀl iÞT iþXiðb iÀm iÞE iþXil i C iÀXi m i C iþXijðk i;jþ1Àk ijÞc i;jÀXil i D iþXi m i D iþXijk ij t ijmij:Since constraints(4b),(5b)and(3)are indepen-dent of each other,RP can be decomposed intotardiness subproblem P T,earliness subproblemP E,and the subproblem to determine f m ij g andf c ij g.P T:minT i P0Xiðw iÀl iÞT i;P E:minE i P0Xiðb iÀm iÞE i;andP MC:minf m ij g;f c ij gL;subject to(3),withLXiðl iÀm iÞC iþXijðk i;jþ1Àk ijÞc i;jÀXiðl iÀm iÞD iþXijk ij t ijmij:Tardiness subproblem P T has an optimal solution:T i¼0if w iÀl i P0,and T i¼þ1if w iÀl i<0.Earliness subproblem P E has an optimal solution:E i¼0if b iÀm i P0,and E i¼þ1if b iÀm i<0.If each operation has only one eligible machinetype to process,m ij will no longer be a decisionvariable and t ijmijwill be a constant.In this case,subproblem P MC can be further decomposed into aset of subproblems,one for each machine type,byregrouping the terms in L according to machinetypes:L¼XhXði;N iÞ2O hðl i(Àm iÞC iþXði;jÞ2O hðk i;jþ1Àk ijÞc ijþXði;jÞ2O hk ij t ijmij)ÀXiðl iÀm iÞD i¼XhXði;jÞ2O h~w ij c ij(þXði;jÞ2O hk ij t ijmij)ÀXiðl iÀm iÞD i¼XhL hÀXiðl iÀm iÞD i;ð8ÞwhereL hXði;jÞ2O h~w ij c ijþXði;jÞ2O hk ij t ijmij;502H.Chen,P.B.Luh/European Journal of Operational Research149(2003)499–512~wij l i Àm i Àk ij ;j ¼N i ;k i ;j þ1Àk ij ;16j <N i :ð9ÞSince each L h in (9)depends only on decision variables related to machine type h ,and Pi ðl i Àm i ÞD i is a constant for a given set of multipliers ff l i g ;f m i g ;f k ij gg ,problem P MC can be decom-posed into a set of subproblems P h ,one for each machine type h ,as follows:P h :min L h ;subject toXijd ijh ðs Þ6M h s ;06s <1:ð10ÞGenerally,some operations may have more thanone eligible machine types to process.In this case,to decompose subproblem P MC into a set of sub-problems,machine types are grouped in a way that any two machine types h 1and h 2with O h 1\O h 2¼;are clustered into the same group,where O h i is the set of operations that can be per-formed on machine type h i as deﬁned in Section 2.Subproblem P MC can then be decomposed into a set of subproblems,each for a group of machine types.If a group of machine types contains only one machine type having only one machine,its corre-sponding subproblem,referred to as the ‘‘machine subproblem’’,is a single machine scheduling prob-lem to minimize the weighted completion time of operations.If the group contains only one machine type but having multiple machines,its corresponding subproblem,referred to as the ‘‘machine type subproblem’’,is an identical par-allel machine scheduling problem [12]to minimizethe same criterion.Otherwise,the group contains more than one machine types.I n this case,its corresponding subproblem,referred to as the ‘‘ma-chine group subproblem’’,is an unrelated parallel machine scheduling problem [12]to minimize the sum of the weighted completion time and the weighted processing time of operations.The sub-problem can be formulated as:P G :min L G X ði ;j Þ2O G~wij c ij þXði ;j Þ2O Gk ij t ijm ij ;subject to Xijd ijh ðs Þ6M h s ;06s <1;h 2G ;ð11Þwhere G is a set of machine types,O G ¼[h 2G O h is the set of operations that can be performed on a machine type belonging to G .The formulations of machine and machine type subproblems can be viewed as two speciﬁc cases of this formulation.For classic job shop scheduling,all subprob-lems of P MC are machine subproblems.Generally,however,the subproblems may contain machine,machine type,and machine group subproblems.The solution of these subproblems will be pre-sented in the next section.3.3.Decomposition and coordination structure As most existing LR approaches,our new ap-proach adopts a two-level decomposition and co-ordination structure as illustrated in Fig.1.AttheH.Chen,P.B.Luh /European Journal of Operational Research 149(2003)499–512503low level,all subproblems are solved for a given set of multipliers,and the set of multipliers are itera-tively adjusted at the high level based on the degree of constraint violation.Simple heuristics are then applied to adjust subproblem solutions to obtain feasible schedules when such an iterative process terminates.The diﬀerences between most existing LR approaches and the new approach are that in our approach,the subproblems are approximately solved by heuristic algorithms,and the recently developed SSG method [20]is used to update Lagrange multipliers at the high level.The reso-lution of the subproblems and the dual problem and the construction of feasible schedules will be presented in the following sections.4.Subproblem solutionputational complexity of the subproblemsFor a given scheduling problem,the subprob-lems of P MC may contain machine subproblems,machine type subproblems,and machine group subproblems.These subproblems have diﬀerent computational complexities.For a machine sub-problem,if the machine is available throughout the scheduling horizon,it is a single machine scheduling problem to minimize the weighted completion time criterion and is polynomially solvable [12].Otherwise,it is NP-complete.Ma-chine type and machine group subproblems,however,are strongly NP-complete even if the number of available machines of each machine type is constant,because identical or unre-lated parallel machine problems to minimize the weighted completion time criterion are NP-complete [12].The NP-complete complexity of machine type subproblems and machine group subproblems implies that optimal algorithms may be not ap-propriate for these subproblems,at least for these with practical sizes [5].Instead,heuristic algo-rithms will be developed in the next subsection.The heuristic resolution of the subproblems is justiﬁed by a recently developed SSG method forthe optimization of Lagrangian dual function [20],as SSG allows approximate optimization of sub-problems under a convergence condition (see the next section for details).The disadvantage of solving subproblems approximately is that the lower bound property of LR is lost,and this is at an exchange of a great reduction in the computa-tion time of our new LR approach.4.2.Solving the subproblemsFor simplicity,all machines are assumed to be available throughout the scheduling horizon in the following description of the heuristic algorithms.For subproblems where some machines are not available during certain time intervals,similar approximate algorithms that take account of ma-chine availability constraints exist after appro-priate modiﬁcations [18].Also for simplicity of statement,two dimensional operation index (i,j )is replaced by one dimensional index i after an ap-propriate index transformation.With the new operation index,the criterion of a subproblem canbe written as J P i ~wi c i þP i k i t im i ,where ~w i ,m i ;c i and t im i correspond to ~wij ,m ij ;c ij and t im ij ,respectively.n the following,three algorithms will be pre-sented for machine,machine type,and machine group subproblems,respectively.Algorithm 1(Algorithm for machine subproblems ).Each machine subproblem is optimally solved by using the Smith Õs weighted shortest processing time (WSPT)rule [12].The rule sequences all op-erations in a non-decreasing order of the ratio t i =~wi ,where t i and ~w i are the processing time and the weight of an operation i ,respectively.The computational complexity of the rule-based algo-rithm is O ðn log ðn ÞÞ,where n is the number of operations.Algorithm 2(Algorithm for machine type subprob-lems ).Each machine type subproblem is approxi-mately solved by using a parameter list scheduling heuristic [12].n the heuristic,all operations are ordered according to a non-decreasing order of the504H.Chen,P.B.Luh /European Journal of Operational Research 149(2003)499–512ratio t i =~w i as in the Smith Õs WSPT rule for machinesubproblems.The heuristic assigns the next oper-ation to the machine that ﬁrst becomes available.In most cases,the WSPT-based heuristic provides a near-optimal schedule for the machine type sub-problem.The computational complexity of the heuristic is O ðmn log ðn ÞÞ,where m is the number of identical machines and n is the number of opera-tions.Algorithm 3(Algorithm for machine group sub-problems ).A machine group subproblem is ap-proximately solved by using a similar parameter list scheduling heuristic.For each machine type h ,all operations that can be processed by it are arranged in a non-decreasing order of the ratio t ih =~wi .This forms H operation lists,one for each machine type.These lists are indexed by machine type,where some operations may appear in more than one list.With these lists,the heuristic cal-culates the cost ~wi c i þk i t ih for each list h and its ﬁrst operation i and selects an operation and a list (a machine type)with the minimum cost as the operation to be dispatched next and the machine type that the operation will be assigned to.I n the calculation,the completion time c i is obtained assuming that operation i will be as-signed to the earliest available machine of type h next.The selected operation is then assigned to the earliest available machine of the selected machine type.After the operation is assigned,its duplications in other lists are removed.This process repeats until all operations are assigned.The computational complexity of the heuristic is O ðmn log ðn ÞÞ,where m is the total number of machines and n is the total number of opera-tions.n few cases,the solution of the relaxed prob-lem obtained by using the above algorithms does not satisfy the convergence condition for SSG at some iterations.I n this case,Algorithm 3is en-hanced by a local search procedure to improve the solution quality.The local search procedure real-locates an operation to another machine type at each step until the condition holds or further im-provement of a solution is impossible.5.Dual problem and the heuristicsn this section,the dual problem is solved by using our recently developed SSG method [20],and feasible schedules are constructed using heu-ristics based on subproblem solutions or the mul-tipliers.Before introducing the method,a propertyof the dual problem is explored to restrict the solution space of the problem.5.1.Property of the dual problemFor an optimal dual solution,all coeﬃcients in the objective functions of P T ;P E ,and P MCð16h 6H Þ,i.e.,w i Àl i ;b i Àm i ,and ~wij should be non-negative for any i and j .Otherwise,by taking T i ¼þ1;E i ¼þ1,or c ij ¼þ1,we have L ¼À1.Since the dual problem is a maximiza-tion problem with a bounded optimal objective value,such multipliers will not be an optimal so-lution.Consequently we have:l i 6w i ;m i 6b i ;ð12a Þk i ;j þ1P k ij ;if 16j <N i ;i ¼1;2;...;N ;ð12b Þl i P m i þk iN i :ð12c ÞLet X be the set of ff l i g ;f m i g ;f k ij gg Õs satisfying constraints (12a)–(12c).We thus have the follow-ing property:Pro erty1.The optimal solution of DP is attained at a point in set X .ntuitively,the multiplier k ij can be interpreted as the price (marginal cost)for one time unit later completion of operation (i ;j À1)or the price for one time unit earlier completion of operation ði ;j Þfor part i .The inequality k i ;j þ1P k ij implies that the late completion of operation ði ;j þ1Þis more costly (crucial)than the late completion of its preceding operation ði ;j Þ.This is because the time slack for the completion of the part before its due date becomes smaller and smaller as time ad-vances.The constraint (12c)is a variant of the constraint (12b)when the ﬁrst and the last oper-ations of a part are considered.The multiplier l iH.Chen,P.B.Luh /European Journal of Operational Research 149(2003)499–512505can be interpreted as the price for one time unit later completion of part i ,while the weight w i is the price for one time unit later delivery of the part.The inequality l i 6w i implies that the late com-pletion of the part is less costly (less crucial)than or as costly (crucial)as the late delivery of the part.This is because one unit later completion of a part does not imply one unit later delivery of the part if the part is completed before its due date.A similar explanation can be given to m i 6b i .5.2.Surrogate subgradient methodLagrangian dual problems for separable integerprogramming problems are commonly solved by using the subgradient method,which requires op-timally solving all subproblems at each iteration to obtain a subgradient direction.This may be time consuming for ones with many subproblems or some hard subproblems.Recently,the surrogate subgradient method has been developed to over-come the diﬃculty [20].n the method,a proper direction can be obtained without optimally solv-ing all the subproblems.In fact,only approxi-mate optimization of one or several subproblems is needed to get a proper ‘‘SSG direction.’’The convergence of the method is proved.Because this method can obtain good directions with much less eﬀort,it is powerful for problems of large size.n the following,a brief introduction of the method is given.Consider a separable integer programming problem described as ðIP Þmin xJ IP ¼X I i ¼1J i ðx i Þ;ð13Þsubject to Ax 6bandx i 2Z n i ;i ¼1;...;I ;ð14Þwhere x ¼ðx 1;x 2;...;x n ÞTis an n Â1decisionvariable with n ¼P Ii ¼1n i and Z is the set of inte-gers.The LR of IP is given byL ðk Þ min x 2ZnXI i ¼1J i ðx i Þ"þk T ðAx Àb Þ#;ð15Þand the Lagrangian dual problem is ðLD Þ:max k P 0L ðk Þ;ð16Þwhere k is a vector of Lagrange multipliers.As an extension of the Lagrangian dual,a sur-rogate dual is introduced:e L ðk ;x Þ XI i ¼1J i ðx i Þ"þk T ðAx Àb Þ#;x 2Z n ;ð17Þand its corresponding SSG is deﬁned as ~gðx Þ Ax Àb :ð18ÞSSG methodStep 0(Initialize).Initialize k 0and minimize all subproblems to obtain x 0i.e.,x 0¼arg min x 2ZnX I i ¼1J i ðx i Þ"þðk 0ÞTðAx Àb Þ#:ð19ÞStep 1(Update multipliers).Given the current point (k k ;x k )at the k th iteration,the Lagrange multipliers are updated according to k k þ1¼k k þs k ~gk ;ð20Þwhere ~gk is the SSG given by ~gk ¼~g ðx k Þ¼Ax k Àb ;ð21Þwith stepsize s k satisfying 0<s k <L Ã Àe L k ~gk 2:ð22ÞHere L Ã¼L ðk Ã)is the optimal objective of dualproblem LD ,e L k ¼e L k ðk k ;x k Þis the surrogate dual at the k th iteration.Step 2(Perform approximate optimization).Given k k þ1perform ‘‘approximate optimization’’to obtain x k þ1such that x k þ1satisﬁes e L ðk k þ1;x k þ1Þ<e L ðk k þ1;x k Þ:ð23ÞIf such an x k þ1cannot be obtained,set x k þ1¼x k .Step 3(Checking stopping criteria).If the cri-teria given by k k k þ1Àk k k <e 1andð24a Þ506H.Chen,P.B.Luh /European Journal of Operational Research 149(2003)499–512。