Bifurcation method for solving multiple positive solutions to boundary value problem of p-Henon

Hopf bifurcationFrom Wikipedia, the free encyclopedia (Redirected from Andronov-Hopf bifurcation )Jump to: navigation , search In the mathematical theory of bifurcations , a Hopf or Poincaré–Andronov–Hopf bifurcation, named after Henri Poincaré, Eberhard Hopf , and Aleksandr Andronov , is a local bifurcation in which a fixed point of a dynamical system loses stability as a pair of complex conjugate eigenvalues of the linearization around the fixed point cross the imaginary axis of the complex plane . Under reasonably generic assumptions about the dynamical system, we can expect to see a small-amplitude limit cycle branching from the fixed point.For a more general survey on Hopf bifurcation and dynamical systems in general, see [1][2][3][4][5].Contents[hide ]● 1 Overview r 1.1 Supercritical / subcritical Hopf bifurcationsr 1.2 Remarks r1.3 Example ● 2 Definition of a Hopf bifurcation ● 3 Routh–Hurwitz criterionr 3.1 Sturm seriesr 3.2 Propositions ● 4 Example● 5 References●6 External links [edit ] Overview[edit ] Supercritical / subcritical Hopf bifurcationsThe limit cycle is orbitally stable if a certain quantity called the first Lyapunov coefficient is negative, and the bifurcation is supercritical. Otherwise it isunstable and the bifurcation is subcritical.The normal form of a Hopf bifurcation is:where z , b are both complex and λ is a parameter. WriteThe number α is called the first Lyapunov coefficient.●If α is negative then there is a stable limit cycle for λ > 0:whereThe bifurcation is then called supercritical.●If α is positive then there is an unstable limit cycle for λ < 0. The bifurcation is called subcritical.[edit ] Remarks The "smallest chemical reaction with Hopf bifurcation" was found in 1995 in Berlin, Germany [6]. The same biochemical system has been used in order to investigate how the existence of a Hopf bifurcation influences our ability to reverse-engineer dynamical systems [7].Under some general hypothesis, in the neighborhood of a Hopf bifurcation, a stable steady point of the system gives birth to a small stable limit cycle . Remark that looking for Hopf bifurcation is not equivalent to looking for stable limit cycles. First, some Hopf bifurcations (e.g. subcritical ones) do not imply the existence of stable limit cycles; second, there may exist limit cycles not related to Hopf bifurcations.[edit ] ExampleThe Hopf bifurcation in the Selkov system(see article). As the parameters change, a limitcycle (in blue) appears out of an unstableequilibrium.Hopf bifurcations occur in the Hodgkin–Huxley model for nerve membrane, the Selkov model of glycolysis , the Belousov–Zhabotinsky reaction , the Lorenz attractor and in the following simpler chemical system called the Brusselator as the parameter B changes:The Selkov model isThe phase portrait illustrating the Hopf bifurcation in the Selkov model is shown on the right. See Strogatz, Steven H. (1994). "Nonlinear Dynamics and Chaos" [1], page 205 for detailed derivation.[edit ] Definition of a Hopf bifurcationThe appearance or the disappearance of a periodic orbit through a local change in the stability properties of a steady point is known as the Hopf bifurcation. The following theorem works with steady points with one pair of conjugate nonzero purely imaginary eigenvalues . It tells the conditions under which this bifurcation phenomenon occurs.Theorem (see section 11.2 of [3]). Let J 0 be the Jacobian of a continuous parametric dynamical system evaluated at a steady point Z eof it. Suppose that all eigenvalues of J 0 have negative real parts except one conjugate nonzero purely imaginary pair. A Hopf bifurcation arises when these two eigenvalues cross the imaginary axis because of a variation of the system parameters.[edit ] Routh–Hurwitz criterionRouth–Hurwitz criterion (section I.13 of [5]) gives necessary conditions so that a Hopf bifurcation occurs. Let us see how one can use concretely this idea [8].[edit ] Sturm series Let be Sturm series associated to a characteristic polynomial P . They can be written in the form:The coefficients c i,0 for i in correspond to what is called Hurwitz determinants [8]. Their definition is related to the associated Hurwitz matrix .[edit ] PropositionsProposition 1. If all the Hurwitz determinants c i ,0 are positive, apart perhaps c k,0 then the associated Jacobian has no pure imaginary eigenvalues.Proposition 2. If all Hurwitz determinants c i ,0 (for all i in are positive, c k " 1,0 = 0 and c k" 2,1 < 0 then all the eigenvalues of the associated Jacobian have negative real parts except a purely imaginary conjugate pair.The conditions that we are looking for so that a Hopf bifurcation occurs (see theorem above) for a parametric continuous dynamical system are given by this last proposition.[edit ] Example Let us consider the classical Van der Pol oscillator written with ordinary differential equations:The Jacobian matrix associated to this system follows:The characteristic polynomial (in λ) of the linearization at (0,0) is equal to:P (λ) = λ2 " μλ + 1.The coefficients are: a 0 = 1,a 1 = " μ,a 2 = 1 The associated Sturm series is:The Sturm polynomials can be written as (here i = 0,1):The above proposition 2 tells that one must have:c 0,0 = 1 > 0,c 1,0 = " μ = 0,c 0,1 = " 1 < 0.Because 1 > 0 and 1 < 0 are obvious, one can conclude that a Hopf bifurcation may occur for Van der Pol oscillator if μ = 0.[edit ] References1. ^ a b Strogatz, Steven H. (1994). Nonlinear Dynamics and Chaos . Addison Wesley publishing company.2. ^ Kuznetsov, Yuri A. (2004). Elements of Applied Bifurcation Theory . New York: Springer-Verlag. ISBN 0-387-21906-4.3. ^ a b Hale, J.; Ko ak, H. (1991). Dynamics and Bifurcations . Texts in Applied Mathematics. 3. New York: Springer-Verlag.4. ^ Guckenheimer, J.; Myers, M.; Sturmfels, B. (1997). "Computing Hopf Bifurcations I". SIAM Journal on Numerical Analysis .5. ^ a b Hairer, E.; Norsett, S. P.; Wanner, G. (1993). Solving ordinary differential equations I: nonstiff problems (Second ed.). New York: Springer-Verlag.6. ^ Wilhelm, T.; Heinrich, R. (1995). "Smallest chemical reaction system with Hopf bifurcation". Journal of Mathematical Chemistry 17 (1): 1–14.doi :10.1007/BF01165134. http://www.fli-leibniz.de/~wilhelm/JMC1995.pdf .7. ^ Kirk, P. D. W.; Toni, T.; Stumpf, MP (2008). "Parameter inference for biochemical systems that undergo a Hopf bifurcation". Biophysical Journal 95 (2):540–549. doi :10.1529/biophysj.107.126086. PMC 2440454. PMID 18456830. /biophysj/pdf/PIIS0006349508702315.pdf .8. ^ a bKahoui, M. E.; Weber, A. (2000). "Deciding Hopf bifurcations by quantifier elimination in a software component architecture". Journal of SymbolicComputation 30 (2): 161–179. doi:10.1006/jsco.1999.0353. [edit] External links● Reaction-diffusion systems● The Hopf Bifurcation● Andronov–Hopf bifurcation page at ScholarpediaCategories: Bifurcation theoryPersonal tools● Log in / create accountNamespaces● Article● DiscussionVariantsViews● Read● Edit● View historyActionsSearchInteractionToolboxPrint/exportLanguages● This page was last modified on 25 May 2011 at 02:56.● Text is available under the Creative Commons Attribution-ShareAlike License; additional terms may apply. See Terms of use for details. Wikipedia is a registered trademark of the W ikimedia Foundation, Inc., a non-profit organization.● Contact us● Privacy policy● About Wikipedia● Disclaimers●●。

解析函数的几种定义及其关系（1）潘传中(达州职业技术学院党办，四川达州635000）摘要：解析函数在复变函数中占有十分重要的地位，但它的定义在不同的著作中形式不一样.给出了解析函数的五种定义，并证明它们彼此等价，从而能更深刻地理解和应用解析函数. 关键词：解析函数；定义；等价Five Definitions and Their Relations of Analytic FunctionPAN Chuan-zhong(Party Office of Dazhou V ocational and Technological College, Dazhou Sichuan 635000, China) Abstract:The Analytic Function occupies a very important position in the Complex Functions, but its definition has various forms in different works. For the sake of getting a better understanding, five kinds of definitions are listed in this article. Their equivalence has also been proved thus to be understood more profoundly and applied more appropriately.Key words: Analytical Function; definition; equivalence关于费马数的一个结果及应用（4）管训贵(泰州师范高等专科学校数理系，江苏泰州225300）摘要：利用费马合数Fn的素因数分解式的一般形式，建立了具有一定应用价值的等式，得到了F5的素因数分解式，并证明了F6是合数.文末还给出了一个猜想.关键词：费马数；合数；素因数分解式；等式；猜想A Result of Fermat Numbers and Its ApplicationGUAN Xun-gui(Mathematics and Polytechnic Department of Taizhou Teachers' College, Taizhou Jiangsu 225300,China)Abstract:In this paper, an equivalence is established by using the generalization of the prime factorization of Fn. It is of great value. We obtain the prime factorization of F5 and F6 is composite. A conjecture is given at end.Key words:Fermat numbers; composite; the prime factorization; equality; conjecture关于人工约束法寻找对偶初始可行解的一个注记（6）董兵1，陈文2(1.中国民航飞行学院计算机学院，四川广汉618307；2.四川文理学院数学与财经系，四川达州635000)摘要：用对偶单纯形法求解线性规划问题，在无法直接求得对偶问题的可行解时，引入人工约束法寻找对偶问题初始可行解.讨论了原问题（LP）与新规划（LPM）解之间的关系，并给出了证明.关键词：线性规划；人工约束法对；偶规划A Note of Artificial Constrained Method Searching for Dual Feasible SolutionDONG Bing, CHEN Wen(puter Department of Civil Aviation Flight University of China, Guanghan Sichuan 618307;2. Mathematics and Finance-Economics Department of SASU, Dazhou Sichuan 635000, China) Abstract:When directly obtaining no feasible solution to the dual problem with the dual simplex method for solving linear programming problems, the artificial constraint method can be introduced to find an initial feasible solution to the dual problem. The relationship between thesolutions of the original problem (LP) and the new planning (LPM) is discussed and the proofs are given in the paper.Key word:Linear programming; artificial constrained method; dual programming二阶常系数线性非齐次微分方程的通解（8）张金战(陇南师范高等专科学校数学系，甘肃成县742500）摘要：在已知二阶常系数齐次微分方程y″+py′+qy=0的一个特解的条件下，讨论了求二阶常系数线性非齐次微分方程y″+py′+qy=f(x)的一个特解的方法，从而根据齐次方程的特征根的不同情形给出了非齐次微分方程的通解公式.关键词：线性微分方程；特解；通解A General Solution to Order 2 Constant Coefficients Non-homogeneous Linear DifferentialEquationZHANG Jin-zhan(Mathematics Department of Longnan Teachers' College, Chengxian Gansu 742500, China) Abstract:On the basis of a special solution to order 2 constant coefficient homogeneous differential equation y″+py′+qy = 0, a special solution is discussed to order 2 constant coefficient non-homogeneous linear differential equation of y″+ py′+ qy = f(x). Following that, the formula of general solution to non-homogeneous differential equation is given.Key word:Linear differential equation; general solution; special solution.关于一例数学建模竞赛题假设合理性的研究（10）王凡彬1，2，陈敏1，唐春梅3(1.内江师范学院数学与信息科学学院，四川内江641112；2. 四川省高等学校数值仿真重点实验室，四川内江6411123.内江师范学院计算机科学系，四川内江641112)摘要：对2001年全国大学生数学建模竞赛C题参考答案的假设中的“奖金在年底发放”提出了异议，认为该假设只是使奖金额达到最大，但方案的可操作性很弱，不符合实际情况.给出了新的符合实际的假设—每年１月１日预留当年奖金，奖金可在一年中任何时候发放，并根据新的假设对Ｃ题重新建立了数学模型，得到了新的投资方案.关键词：数学建模；竞赛；假设；线性规划Assumption Rationality of a Case in 2001 National Undergraduate Math Modeling ContestWANG Fan-bin1, 2, CHEN Min1, TANG Chun-mei3(1. Mathematics and Information Science College, Neijiang Normal University, Neijiang Sichuan641112;2. Key Laboratory of Numerical Simulation of Sichuan Higher Education, Neijiang Sichuan641112, China)Abstract:The assumption of the reference answer to Question C is questioned that “The bonus is paid at the e nd of the year” in 2001 National Undergraduate Math Modeling Contest. The writer of this paper argues that the assumption just only maximizes the prize money, but its operation is very weak and inconsistent with the actual situation. Realistic assumptions is given: the money is set aside on January 1and the bonus can be paid at any time of the year. According to the new assumption, a new math model of Question C is established and a new investment program is achieved.Key words:Math Modeling Contest; contest; assumption; linear programming具有Leslie-Gower反应的离散捕食系统的稳定性和分支分析（13）张莉敏(四川文理学院数学与财经系，四川达州635000）摘要：研究了一类用向前欧拉法获得的具有Leslie-Gower反应类型的离散捕食系统的动力学行为.利用Jury判据，探讨了系统的渐进稳定性；利用分支理论和中心流型定理，证明了系统在一定条件下存在flip分支.关键词：向前欧拉法；离散捕食-食饵系统；flip分支Stability and Bifurcation in a Discrete Predator-prey System with Leslie-Gower TypeZHANG Li-min(Mathematics and Finance-Economics Department of SASU, Dazhou Sichuan 635000, China) Abstract:The dynamic behavior of a discrete predator-prey system obtained by forward Euler method with Leslie-Gower type is investigated. The Criterion Jury is adopted to analyze the asymptotical stability. The center manifold theory and bifurcation theorem can prove that the flip bifurcation exists in a certain condition.Key words:Forward Euler method; discrete predator-prey system; flip bifurcation一类捕食者—食饵模型的周期解与稳定性（16）马丽蓉(四川民族学院数学系，四川康定626001）摘要：研究了一类含扩散与时滞捕食者—食饵模型，利用上下解及比较原理，证明了在一定条件下该模型的零平衡态及半平凡周期解的全局稳定性，并获得了这个系统具有一对周期拟解的充分条件，且对任意的非负初值函数，这对周期拟解构成的区间是此系统的一个吸引子. 关键词：扩散；时滞；捕食者-食饵模型；上下解；周期解Stability and Periodic Solution for a Prey-Predator ModelMA Li-rong(Mathematics Department of Sichuan Nationalities College, Kangding Sichuan 626001, China) Abstract:The existence and stability of periodic solution in Prey-Predator model with diffusion，time-delay are investigated by constructing a pair of upper and lower solutions and comparison principle. It is shown that under some appropriate conditions the trivial solution and semi-trivial periodic solution of the model are globally asymptotically stable, the models have a pair of periodic quasi-solutions and the sector between the quasi-solutions is an attractor of the model with respect to every nonnegative initial function.Key words:Diffusion; time delays; prey-Predator model; upper and lower solutions; periodic solution完全单半群上同余的另一刻画（20）罗肖强(四川文理学院数学与财经系，四川达州635000）摘要：在完全单半群的幂等元集E（S）上构建正规等价τ与核子集K，证明了K与τ的相容性，以及正规子群之间的蕴含关系，迹类与正规子群格之间的同构，核类与幂等元集E(S)上的正规等价格同构.关键词：完全单半群；正规等价；迹核同余；相容性Another Depiction of Congruence on Completely Simple Semi-groupsLUO Xiao-qiang(Mathematics and Finance-Economics Department of SASU, Dazhou Sichuan 635000, China) Abstract:On the idempotent element set E (S) of the completely simple semi-groups, the normalequivalence τ and the kernel subset K are structured, and it is showed that the K and τ is compatible, their normal subgroups are implicative, the trace classes and normal subgroup lattices are isomorphic, and the kernel classes and the normal equivalence lattices within the idempotent element set are isomorphic, too.Key words:Completely simple semi-groups; normal equivalence; the trace-kernel congruence; compatible property利率具有二阶自回归相依结构的破产问题（23）李爱民1，涂庆伟2(1.四川文理学院数学与财经系，四川达州635000；2.江苏工业学院数理学院，江苏常州213016 )摘要：研究了利率具有二阶自回归相依结构的风险模型，通过递推关系，得到了破产前盈余分布和首达某一水平x的时间分布的积分方程.关键词：二阶自回归相依模型；破产前盈余分布；递推公式Bankruptcy Problems under Interest Rates with Order 2 Autoregressive StructureLI Ai-min, TU Qing-wei( 1. Mathematics and Finance-Economics Department of SASU, Dazhou Sichuan 635000, China；2. Mathematics Department of Jiangsu Polytechnic University, Changzhou Jiangsu 213016,China)Abstract:The risk model is studied that the interest rate is dependent upon the second autoregressive structure. By recursive relations, the integrate equations are derived of the surplus distribution before bankruptcy and the time X distribution that reaches a given level for the first time.Key words:Order 2 autoregressive model; surplus distribution before bankruptcy; recursive formula科技期刊影响因子的偏差分析（26）周兴旺(四川大学学报编辑部，四川成都610064）摘要：作为国际公认的衡量期刊影响力以及能否入选SCI数据库的一个重要指标，影响因子备受科研管理者、科研工作者和办刊者的重视，而影响因子又容易出现偏差。

Mooney-Rivilin压缩杆模型的精确解及其动力学性质王婧【摘要】利用动力系统理论中的积分分支法研究Mooney-Rivilin压缩杆模型的行波解,获得了包括尖孤子解、爆破波解、周期波解、亮孤子解和暗孤子解等各种精确行波解.进一步讨论了这些精确行波解随时间演化的动力学行为,并利用Maple 软件绘出了具有代表性的精确行波解随时间演化的坐标图形.与现有文献的结果相比,本文所获得的精确解比较新颖.【期刊名称】《玉溪师范学院学报》【年(卷),期】2017(033)012【总页数】6页(P15-20)【关键词】积分分支法;精确行进波;孤立波解;尖孤子解【作者】王婧【作者单位】重庆师范大学数学科学学院,重庆 401331【正文语种】中文【中图分类】TQ330自然科学中的许多非线性问题，都可以用非线性偏微分方程和发展方程来建模,且非线性偏微分方程和发展方程的解可以用来解释模型的某些物理现象,特别是质点随时间演化的各种动力学行为和动力学性质.有时为了能够更加精准地了解模型所代表的事物内在变化的规律，我们就要去求其相应模型的精确解.目前,在精确求解非线性偏微分方程和发展方程方面已获得的一些有效的方法，如：反散射变换[1,2]、达布变换[3,4]、双线性方法和多元线性方法[5,6]、齐次平衡法[7～9]、Lie 群法[10,11]、平面动力系统分支理论[12]等.在文献[13,14]中，Rui等在平面动力系统分支理论[12]的基础上，提出了一个比较简化的方法叫做积分分支方法，该方法与平面动力系统分支理论中方法相比，绕开了相对复杂的平面像图分析过程，采取因式分解技巧与直接积分相结合的方法来获得非线性偏微分方程和发展方程的各种精确行波解并进一步讨论解的各种动力学行为和动力学性质，以此来解释模型所代表的问题中各种非线性物理现象.本文将利用积分分支法与因式分解相结合的方法来调查下列非线性Mooney-Rivilin压缩杆模型ut-uxxt+3uux-λ(2uxuxx+uuxxx)+βux=0(1)的各种精确解，并进一步讨论它们的动力学行为和动力学性质.方程(1)是戴辉辉从可压缩的超弹性材料中导出的弹性杆模型，又称Mooney-Rivilin压缩杆模型，其中u=u(x,t),而参数λ为材料参数依赖于预应力，γ为一常数[15].特别地，当λ=1时，方程(1)就可约化成著名的潜水波模型Camassa-Holm方程，曾经被许多研究人员研究过.同样，这个可压缩的超弹性杆模型(1)也被很多研究人员研究过它的行波的间断和爆破现象，详细报道见文献[16]中所引用的大量文献.2 模型的精确行波解及动力学性质为了获得方程(1)的各种精确行波解，我们做下列行波变换u(x,t)=φ(ξ),ξ=lx-ct(2)其中l为波数，c为波速度，二者均为非零实常数.将(2)式代入(1)式后整理得：(βl-c)φ′+cl2φ‴+3lφφ′-2λl3φ′φ″-λl3φφ‴=0,(3)其中φ′=dφ/dξ将(3)式积分一次得：g+(βl-c)φ-0.5λl3(φ′)2+(cl2-λl3φ)φ″+1.5lφ2=0,(4)其中g为积分常数.当φ≠c/(λl)时，方程(4)化为φ″=[g+(βl-c)φ-0.5λl3(φ′)2+1.5lφ2]/(λl3φ-cl2).(5)令φ′≡dφ/dξ=y,则方程(5)可化为下列平面动力系统：dφ/dξ=y,dy/dξ=[g+(βl-c)φ-0.5λl3y2+1.5lφ2]/(λl3φ-cl2),(6)系统(6)是一个奇异系统，因为当φ=c/(λl)时，dy/dξ无定义，由此我们称φ=c/(λl)为系统(6)的奇异直线.显然当φ=c/(λl)时，系统(6)与方程(4)不等价，然而，φ=c/(λl)也是方程(4)的一个解(一个常数解).为了获得与方程(4)完全等价的平面动力系统，我们做下列尺度变换dξ=(λl3φ-cl2)dτ,(7)其中τ为参数.则系统(6)在变换(7)下被化成一个规则的二维平面动力系统：=y(λl3φφ-0.5λl3y2+1.5lφ2.(8)显然无论变量φ如何变化，系统(8)始终与方程(4)等价.通过简单计算后不难发现系统(6)和(8)具有相同的首次积分，即0.5lφ3+0.5(βl-c)φ2-0.5λl3y2φ+0.5cl2y2+gφ=h,(9)其中h为新的积分常数.记函数H(φ,y)=0.5lφ3+0.5(βl-c)φ2-0.5λl3y2φ+0.5cl2y2+gφ=h.(10)容易验证函数H(φ,y)满足∂H/∂φ=dy/dτ,∂H/∂y=-dφ/dτ所以系统(8)是一个哈密顿系统，满足哈密顿守恒律.显然O(0,0)原点始终是系统(8)的其中一个平衡点.利用(10)式，不难计算出原点的Hamilton量为h0=H(0,0)=0.根据积分常数和哈密顿量是否为零的情况，下面讨论方程(1)的精确行波解.首先，由方程(9)可解得(11)情形1 当两个积分常数g=0,h=0时，将(11)式代入系统(8)的第一式后化成=±dτ(12)在计算的过程中，为方便起见，记(12)式等号右边的正负号为ε,那么在不同的参数条件下积分(12)后换回变量u，并将积分所得的结果分别代入变换式(7)中再次积分后联立u和ξ两式，我们便得到下列方程(1)的12种参数形式的精确行波解：(13)(14)(15)(16)(17)(18)(19)(20)(21)(22)(23)(24)在以上这些参数形式的精确行波解中，解(13)(14)(15)为孤立波解，解(17)(18)(20)为周期波解，解(23)和(24)为有理解.为了直观地显示以上解的动力学行为，以其中部分解为例画出其坐标图，见图1.图1 参数条件g=0,h=0下各种参数型精确行波解的坐标图情形2 当g≠0,h=0,λ>0时，将(11)式代入(7)式和(8)式的第一个式子可得φ(25)其中φ3=c/(λl),φ4=0.当积分常数满足g=(c-βl)2/(8l)条件且c=λβl/(λ-2)时，则有φ1,2=c/(λl),那么(25)式可化成(26)完成(26)式的积分后换回变量，我们获得方程(1)的一个显函数形式的精确行波解：(27)当积分常数满足g=(c-βl)2/8l条件，但c≠λβl/(λ-2)时，则仅有φ1=φ2=(c-βl)/2l,那么(26)可化为φ(28)完成(30)式的积分后换回变量，我们获得方程(1)的一个隐函数形式的精确行波解：(29)为了直观地显示精确行波解(27)和(29)的动力学行为，我们画出了它们的坐标图,如图2所示.情形3 当积分常数g≠0,h≠0时，将(11)式直接代入(8)式中的第一个式子化简可得(30)(a)在特定的参数条件下(30)式可分解为φ(31)图2 由解(27)和(29)定义的U型波当φ1>φ2>φ3≥0>φ,完成(31)式积分后换回相应变量，可获得方程(1)的一种隐函数形式的解：其中为第三类椭圆积分函数且式中其他参数为(b)在特定的参数条件下，(30)又可以分解为φ(33)其中φ1>φ2≥φ>φ3>0.类似地，完成(33)式的积分后换回相应变量，我们得到方程(1)的另一个隐函数形式的精确行波解：(34)其中为第三类椭圆积分函数且式中其他参数为该文由芮伟国教授指导完成.参考文献：[1]V.E.Zakharov and A.B.Shabat,Exact theory of two-dimensional self-focusing and one dimensional self-modulation of waves in nonlinear media[J].JETP,1972(34):62-69.[2]C.S.Garder,J.M.Greene,M.D.Kruskal,R.M.Mirura,Method for solving the Korteweg-de -Vries equation[J].Physical Review Letters,1967(19):1095-1097.[3]C.S.Garder,J.M.Greene,M.D.Kruskal,R.M.Mirura,Korteweg-de Vries equation and generalizations[J].VI.Methods for exactsolution,Communications on Pure and Applied Mathematics,1974,27(1):97-133.[4]V.B.Matveev,M.A.Salle,Darboux transformations andSolitons[M].Berlin:Springer,1991.[5]I.T.Habibullin, Backlund transformation and integrable initial-boundary value problems[J].Mathematical notes of the Academy of Sciences of the USSR,1991,49(4):418-423.[6]R.Hirota,Exact solution of the Korteweg-de Vries equation for multiple collisions of solitons[J].Physical Review Letters,1971(27):1192-1194. [7]A.M.Wazwaz,The tanh method for traveling wave solutions of nonlinear equations[J].Applied Mathematics and Computation,2004,154(3):713-723.[8]M.Wang,Y.Zhou,Z.Li,Application of a homogeneous balance method to exact solutions of nonlinear equations in mathematical physics[J].Physics Letters A,1996,216(1):67-75.[9]E.Fan,H.Zhang,A note on the homogeneous balance method[J].Physics Letters A,1998,246(5):403-406.[10]E.D.Belokolos,A.I.Bobenko,Algebro-geometric approach to nonlinear integrable Equations[M].Berlin:Springer-Verlag,1994.[11]N.C.Freeman,J.J.C.Nimmo,Soliton solutions of the Korteweg de Vries and the Kadomtsev-Petviashvili equations: the Wronskiantechnique[J].Proceedings of the Royal Society of London A: Mathematical, Physical and Engineering Sciences.The Royal Society,1977(389):319-329.[12]J.Li,Z.Liu,Smooth and non-smooth travelling waves in nonlinearly dispersive equation[J].Applied Mathematical Modeling,2000(25):41-56. [13]W.Ru,B.He,Y.Long,et al.,The integral bifurcation method and its application for solving a family of third-order dispersive PDEs[J].Nonlinear Analysis: Theory, Methods and Applications,2008,69(4):1256-1267.[14]W.Rui,B.He,S.Xie,et al.Application of the integral bifurcation method for solving modified Camassa-Holm and Degasperis-Procesiequations[J].Nonlinear Analysis:Theory, Methods andApplications,2009,71(7):3459-3470.[15]H.H.Dai,Model equations for nonlinear dispersive waves in a compressible Mooney-Rivlin Rod[J].Acta Mech.1998(127):193-207. [16]F.Ma,Y.Liu,C.Qu,Wave-breaking phenomena for the nonlocal Whitham-type equations[J].Journal of Differential Equations,2016(261):6029-6054.。

P roce dia Compute r Scie nce 22 ( 2013 )512 – 5201877-0509 © 2013 The Authors. Published by Elsevier B.V.Selection and peer-review under responsibility of KES International doi: 10.1016/j.procs.2013.09.130ScienceDirectAvailable online at 513M asaru Noda / P rocedia Computer Science 22 ( 2013 )512 – 520NomenclatureBD blockage degree [%]rate [m3/s]f flowI set of microchannel [-]number [-]i channell length [m]n number [-]vector [Pa]dataP pressurep pressure [Pa]coefficient [-]R correlationS set of sensor locations [-]location [-]s sensorvelocity [m/s]v inletw width [m]z depth [m]ȝviscosity [Pa䡡s]<subscripts>conditionB blockedbot bottomC channelwall)(ChannelF finM manifoldconditionN normalO outlettop topMicrochannels are prone to blockage due to side reactions or contamination from raw materials when they are operated for a long period. Blockage in microchannels causes poor uniformity in the residence time distribution among them, leading to degraded product quality. Blockage in the microchannels of microreactors is a serious problem that limits their practicality. It is therefore essential to detect and identify blockage locations to ensure more effective and stable microreactor operation.Kano et al. [3] proposed data-based and model-based blockage diagnosis methods using temperature sensors that identify a blockage location in stacked microchemical processes. The data-based method compares the ratios of temperature differences between normal and abnormal operating conditions at one sensor to those at the other sensor. The simulation results showed that this method could diagnose the blockage location successfully. However, it might not work if the blockage does not affect the temperature in a microreactor due to its high surface/volume ratio.Tanaka et al. [6] developed a blockage detection and diagnosis system for parallelized microreactors with split-and-recombine-type flow distributors. This system can isolate a blocked microreactor with a small number of flow sensors. Yamamoto et al. [7] proposed a method that uses pressure sensors instead of temperature sensors, in which the blockage location is identified by comparing measured pressure distribution data with prepared pressure distribution data calculated by computational fluid dynamic (CFD) simulation when a blockage occurs. Simulation results showed that these methods could diagnose a single blockage location using514M asaru Noda / P rocedia Computer Science 22 ( 2013 )512 – 520515 M asaru Noda / P rocedia Computer Science 22 ( 2013 )512 – 520516M asaru Noda / P rocedia Computer Science 22 ( 2013 )512 – 520517M asaru Noda / P rocedia Computer Science 22 ( 2013 ) 512 – 520Table 1. Geometric parameters and operation conditions of microreactorName Parameters Value UnitNumber of channels n C 10 - Channel width w C 100 ȝm Channel depth z C 100 ȝm Channel length l C 20 mm Width of fin w F 284 ȝm Viscosity ȝ 0.1 Pa·sWidth of manifoldw M,top 1.0 mmw M,bot 5.0 mmInlet velocity v 0.01 m/s Outlet pressure p O 101.3 kPaTwo blockage diagnosis problems, Case 1 and Case 2, were considered to assess the proposed method’s effectiveness, where four pressure sensor configurations (Conf. A–D) were implemented. In case studies, one of the blocked channels, i 1, was fixed at one or two to reduce the total number of CFD simulations.Conf. A: n S = 10, S = {1, 2, 3, 4, 5, 6, 7, 8, 9, 10} Conf. B: n S = 6, S = {1, 3, 5, 6, 8, 10} Conf. C: n S = 5, S = {1, 3, 5, 7, 9} Conf. D: n S = 4, S = {1, 4, 7, 10}Case 1: ),,,1(%35%,35%,153********i i i I i i i BD BD BD i i iCase 2: ),,,2(%35%,35%,153********i i i I i i i BD BD BD i i i4.2. Diagnosis resultsTen CFD simulations were executed to obtain pressure distribution data under 50% blockage in eachmicrochannel. The simulation results were then used as a data set in a database. Then, sixty-four CFD simulations were executed to obtain a validation data set for thirty-six combinations of three blocked channels for Case 1 and twenty-eight combinations of three blocked channels for Case 2.Simulation results for Cases 1 and 2 are respectively summarized in Tables 2 and 3, in which the combinations of identified numbers of three blocked channels are shown. Several misdiagnoses occurred, which are denoted with underlines in the tables. Table 4 lists the percentages of correct diagnoses for sensor configurations A-D. The method perfectly identified the three blockage locations for both Cases 1 and 2 when using ten pressure sensors (Conf. A). The performance of the proposed method deteriorated when the number of pressure sensors decreased. However, when using six or five sensors (Conf. B and Conf. C), the method still correctly identified the three blocked microchannels in most cases: 100 and 94.4% for Case 1 and 96.4 and 89.3% Case 2. The average correct diagnosis percentage was 36.1% for Case 1 and 21.4% for Case 2 when the number of sensors was four. This is because too small a number of sensors make it impossible to distinguish the changes in patterns of pressure distribution with regard to blockage locations.The average correct diagnosis percentages using six sensors (Conf. B) was 94.5% when the number of blocked channels are two [4]. These results demonstrate that a small database with pressure distribution data for the blockage in a single microchannel can correctly identify multiple blockage locations in a microreactor when the number of pressure sensors is more than half that of microchannels.518M asaru Noda / P rocedia Computer Science 22 ( 2013 )512 – 520 Table 2. Diagnosis results of case 1No. BlockedchannelsSensor configurations (n S)A(10)B(6) C(5) D(4)1 1, 2, 3 1, 2, 31, 2, 31, 2, 32, 3, 82 1, 2, 4 1, 2, 41, 2, 41, 2, 41, 4, 53 1, 2, 5 1, 2, 51, 2, 51, 2, 51, 2, 74 1, 2, 6 1, 2, 61, 2, 61, 2, 61, 5, 65 1, 2, 7 1, 2, 71, 2, 71, 2, 71, 2, 76 1, 2, 8 1, 2, 81, 2, 81, 2, 82, 3, 87 1, 2, 9 1, 2, 91, 2, 91, 2, 92, 3, 98 1, 2, 10 1, 2, 101, 2, 101, 2, 101, 2, 109 1, 3, 4 1, 3, 41, 3, 41, 3, 42, 4, 810 1, 3, 5 1, 3, 51, 3, 51, 3, 52, 5, 811 1, 3, 6 1, 3, 61, 3, 62, 3, 61, 5, 712 1, 3, 7 1, 3, 71, 3, 71, 3, 71, 9, 1013 1, 3, 8 1, 3, 81, 3, 81, 3, 82, 6, 814 1, 3, 9 1, 3, 91, 3, 91, 3, 92, 9, 1015 1, 3, 10 1, 3, 101, 3, 101, 3, 101, 3, 1016 1, 4, 5 1, 4, 51, 4, 51, 4, 52, 3, 817 1, 4, 6 1, 4, 61, 4, 61, 4, 61, 3, 1018 1, 4, 7 1, 4, 71, 4, 71, 4, 71, 4, 719 1, 4, 8 1, 4, 81, 4, 81, 4, 86, 7, 1020 1, 4, 9 1, 4, 91, 4, 91, 4, 96, 7, 921 1, 4, 10 1, 4, 101, 4, 101, 4, 101, 4, 1022 1, 5, 6 1, 5, 61, 5, 63, 4, 82, 6, 823 1, 5, 7 1, 5, 71, 5, 71, 5, 71, 5, 724 1, 5, 8 1, 5, 81, 5, 81, 5, 82, 9, 1025 1, 5, 9 1, 5, 91, 5, 91, 5, 92, 6, 926 1, 5, 10 1, 5, 101, 5, 101, 5, 102, 8, 927 1, 6, 7 1, 6, 71, 6, 71, 6, 73, 4, 728 1, 6, 8 1, 6, 81, 6, 81, 6, 81, 6, 829 1, 6, 9 1, 6, 91, 6, 91, 6, 91, 6, 930 1, 6, 10 1, 6, 101, 6, 101, 6, 101, 6, 1031 1, 7, 8 1, 7, 81, 7, 81, 7, 81, 7, 832 1, 7, 9 1, 7, 91, 7, 91, 7, 93, 4, 533 1, 7, 10 1, 7, 101, 7, 101, 7, 101, 6, 834 1, 8, 9 1, 8, 91, 8, 91, 8, 91, 8, 935 1, 8, 10 1, 8, 101, 8, 101, 8, 101, 8, 1036 1, 9, 10 1, 9, 101, 9, 101, 9, 101, 9, 10519M asaru Noda / P rocedia Computer Science 22 ( 2013 )512 – 520Table 3. Diagnosis results of case 2No. BlockedchannelsSensor configurations (n S)A(10)B(6) C(5) D(4)1 2, 3, 4 2, 3, 41, 3, 42, 3, 42, 4, 52 2, 3, 5 2, 3, 52, 3, 52, 3, 52, 4, 53 2, 3, 6 2, 3, 62, 3, 62, 3, 61, 3, 64 2, 3, 7 2, 3, 72, 3, 72, 3, 72, 3, 75 2, 3, 8 2, 3, 82, 3, 82, 3, 83, 4, 86 2, 3, 9 2, 3, 92, 3, 92, 3, 92, 3, 97 2, 3, 10 2, 3, 102, 3, 102, 3, 101, 3, 78 2, 4, 5 2, 4, 52, 4, 52, 4, 53, 4, 69 2, 4, 6 2, 4, 62, 4, 62, 4, 63, 4, 610 2, 4, 7 2, 4, 72, 4, 72, 4, 72, 4, 711 2, 4, 8 2, 4, 82, 4, 82, 4, 84, 6, 712 2, 4, 9 2, 4, 92, 4, 92, 4, 92, 3, 1013 2, 4, 10 2, 4, 102, 4, 102, 4, 102, 4, 1014 2, 5, 6 2, 5, 62, 5, 62, 5, 61, 3, 515 2, 5, 7 2, 5, 72, 5, 72, 5, 74, 6, 916 2, 5, 8 2, 5, 82, 5, 82, 5, 81, 2, 417 2, 5, 9 2, 5, 92, 5, 92, 5, 97, 9, 1018 2, 5, 10 2, 5, 102, 5, 102, 5, 102, 5, 1019 2, 6, 7 2, 6, 72, 6, 72, 6, 76, 9, 1020 2, 6, 8 2, 6, 82, 6, 82, 6, 84, 6, 721 2, 6, 9 2, 6, 92, 6, 92, 6, 94, 5, 1022 2, 6, 10 2, 6, 102, 6, 102, 6, 104, 7, 823 2, 7, 8 2, 7, 82, 7, 82, 7, 84, 6, 724 2, 7, 9 2, 7, 92, 7, 92, 7, 96, 9, 1025 2, 7, 10 2, 7, 102, 7, 102, 7, 102, 7, 1026 2, 8, 9 2, 8, 92, 8, 92, 5, 101, 7, 827 2, 8, 10 2, 8, 102, 8, 102, 8, 92, 6, 928 2, 9, 10 2, 9, 102, 9, 102, 8, 91, 8, 9Table 4. Diagnosis results of proposed methodConf.(n S) A(10) B(6) C(5) D(4)Case 1 100 % 100 % 94.4 % 36.1 %Case 2 100 % 96.4 % 89.3% 21.4%5.ConclusionA data-based identification method that can identify blockages in multiple microchannels was proposed. To prevent a combinatorial explosion of CFD simulations for database construction, the proposed method identifies blocked locations by using only pressure distribution data when a single channel is blocked. The520M asaru Noda / P rocedia Computer Science 22 ( 2013 )512 – 520results of CFD simulations showed that the method could accurately identify three blockage locations using fewer pressure sensors than there were microchannels. We focused on a case in which three microchannels were blocked, but our method can be easily applied to cases in which there are four or more blocked microchannels without incurring any combinatorial explosion of CFD simulations for database construction. In the proposed method, it is assumed that the number of blocked channels is known when identifying the locations of blocked channels. I intend to work on developing a method to estimate the number of blocked channels in my future research.References[1] Ehrfeld, W., V. Hessel, and H. Lowe, 2000, Microreactors: New Technology for Modern Chemistry, Wiley, VCH, Weinheim.[2] Hessel, V. H., H. Lowe, A. Muller, and G. Kolb, 2005, Chemical Micro Process Engineering, Viley-VCH Verlag, Weinheim.[3] Kano, M., T. Fujioka, O. Tonomura, S. Hasebe, and M. Noda, 2007, Data-Based and Model-Based Blockage Diagnosis for StackedMicrochemical Processes, Chem. Eng. Sci., 62, 1073-1080[4] Noda, M., N. Sakamoto, 2012a, Blockage Diagnosis Method of Multiple Blocked Channels in a Microreactor, Journal of ChemicalEngineering of Japan, 45, 228-232[5] Noda, M., N. Sakamoto, 2012b, Estimation Method of Blockage Degrees of Multiple Channels in a Microreactor, Journal ofChecamical Engineering of Japan, 45, 498-503[6] Tanaka, Y., O. Tonomura, K. Isozaki, and S. Hasebe, 2011, Detection and Diagnosis of Blockage in Parallelized Microreactors, ChemEng. J., 167, 483-489[7] Yamamoto, R., M. Noda, and H. Nishitani, 2009, Blockage Diagnosis of Microreactors by Using Pressure Sensors, Kagaku KogakuRonbunshu, 35(1), 170–176 (in Japanese)。

基于虚拟现实的血管内介入手术三维导丝运动模拟周正东;Pascal Haigron;Vincent Guilloux;Antoine Lucas【摘要】导管和导丝在血管中的运动模拟在介入手术训练、计划及术中辅助治疗中具有重要意义.本文提出了一种快速有效的碰撞消除方法,开发了实时三维介入手术模拟系统,以模拟导管或导丝在实际血管中的运动行为.采用OpenGL图形库检测导管或导丝与血管壁之间的碰撞,通过几何分析和旋转角传播方法消除碰撞,最后对导管或导丝模型施加松弛过程,使其状态与实际状态更加吻合.实验结果表明,导管或导丝模型的运动状态与给定的材料参数密切相关,松弛过程使其状态更加自然,模拟可满足实时要求,方法可靠有效.【期刊名称】《南京航空航天大学学报（英文版）》【年(卷),期】2010(027)001【总页数】8页(P62-69)【关键词】导管;虚拟现实;导丝;多体模型;血管内介入【作者】周正东;Pascal Haigron;Vincent Guilloux;Antoine Lucas【作者单位】南京航空航天大学材料科学与技术学院,南京,210016,中国;法国国家健康与医学研究院,雷恩,35042,法国;雷恩第一大学信号与图像处理实验室,雷恩,35042,法国;中法生物医学信息研究中心,雷思,35042,法国;法国国家健康与医学研究院,雷恩,35042,法国;雷恩第一大学信号与图像处理实验室,雷恩,35042,法国;法国国家健康与医学研究院,雷恩,35042,法国;雷恩第一大学信号与图像处理实验室,雷恩,35042,法国;法国国家健康与医学研究院临床高新技术研究中心,雷恩第一大学附属医院,雷恩,35033,法国【正文语种】中文【中图分类】TP391.41;R445.39INTRODUCTIONIn the past decades,the rapid advances of minimally invasive surgery have led to its wide usage in clinic due to many advantages,such as fast recovery and short stay at hospital.Many devices are involved in such surgery,where guide wires and catheters are the most oftenused.However,the technique is complicated and physicians require extensive training periods to achieve the competency.During the intervention,guide wires and catheters are advanced throughpushing,pulling and twisting by the physician to reach the target location (aneurysm and stenosis).The type selection of the guide wire and the catheter as well as the navigation of them to the target location are all-important for a successful endovascular intervention. The specific patient vascular requires the guide wire and the catheter with specific properties,such as shape,strength,torque,and elasticity,etc.Currently,the selection of the guide wire and the catheter is a difficult task requiring strong clinical expertise.The virtual reality simulation of the intervention can provide a training environment and be useful for the selection of the guide wire and the catheter.By defining the material properties of the guide wire and the catheter,the behavior of the guide wire or the catheter motion can be realistically simulated in a specific patient arterynetwork,thus helping the physician improve the selection of the suitable guide wire and catheter,and make good operation planning. The challenges of the simulation include the physical realistic modeling of the guide wire,the catheter and the vascular;the calculation of the guide wire or catheter motion;tradeoffs among physical and visual reality;the computation cost;and the demands of real-time interaction.The techniques of the endovascular intervention simulation have been investigated by several groups[1-6].The algorithms can be classified as physical or geometrical methods. Geometrical methods,such as splines and snakes,are based on a simplified physical principle to achieve the reality-like results. They are fast yet without physical properties.The main physical approaches to soft tissue modeling are the mass-spring,multi-body dynamics and the finite element modeling(FEM)methods.FEM is the most realistic method for modeling the tissue deformable behavior if the properties of the model are correctly chosen.It describes a shape as a set of basic geometrical elements and the model is defined by thechoice of its elements,its shape function,and other global parameters[3].FEM treats a problem in acontinuous manner,but solves the problem for each element in a discreteway.It requires enormous computation,thus being hard for the real-time simulation.Ref.[4]proposed a real-time deformation methodology of catheters and guide wires during the navigation inside the vascular by combining a real-time incremental finite element model,an optimization strategy based on the substructure decomposition,and a new method for handling thecollision.Ref.[5]proposed an FEM modeling of the guide wire and the vascular structure to simulate the tool-vessel interaction on patientspecific vascular models extracted in near realtime,and provided a preoperative knowledge of the navigation and the behavior of instruments inside the specific patient vascular for planning applications.The mass spring is a common method for real-time simulations,where masses are assigned to vertices(nodes)and a set of springs are allocated to connect vertices.Mass spring methods are easy to be built.Although the simulation levels are not as high as those for FEM,they can produce acceptable real-timesimulations.However,these models are not suitable for modeling catheters and guide wires because the constraint that the length of theguide wire remains constant is not easily incorporated in the deformation modeling[1].The multi-body dynamics is a suitable technique for solving the simulation problem[6].It is most frequently used in the robot control area,when the robot consists of several parts connected by joints.The Featherstone algorithm[7]is used to calculate the acceleration of these parts resulting from forces.However,a guide wireis a flexible device.The calculation of guide wire movements involves the flexibility of the guide wire and the spring energy associated with bending.Thus,the Featherstone algorithm,which does not incorporate the spring energy,is not suitable for the guide wire and the catheter simulation.Although many problems for the guide wire and catheter simulation are addressed in previous work,challenges still exist,especially for seeking ahigher level of the fidelity and the accuracy in the simulation while maintaining the real-time computational performance[1].Based on Ref.[8],this paper focuses on the three-dimensional(3-D)real-time simulation of minimally invasive endovascular intervention to provide the pre-operative knowledge of the guide wire and catheter behavior inside a specific vascular for the surgeon,thus being helpful to choose the guide wire and the catheter with suitable properties for specific patient data.The vascular is segmented from computed tomography(CT)data and represented by the mesh surface.The guide wire/catheter is modeled as a multi-body,and the properties are defined by its intrinsic characteristics with the strength and the elasticity.The motion of theguide wireand the catheter inside the vascular is guided by the collision detection and the cancellation algorithm. The scheme of the navigation procedure is shown in Fig.1.Open graphics library(OpenGL)orthographic camera is used for the collision detection between the guide wire/catheter and the vascular wall,while ageometry analysis method is used to find the right motion direction of the guide wire/catheter.Finally,a relaxation procedure is applied to the model to achieve more realistic status.Fig.1 Scheme of navigation procedure1 METHOD AND MATERIALS1.1 Patient data descriptionThe precise knowledge of the geometrical parameters for arteries and lesions is required for the endovascular intervention simulation. The 3-D geometrical description of the vascular inner surface is assumed to be rigidbased on a virtual angioscopy likeprocess applied to CTdata.In the virtual exploratory navigation framework, the virtual angioscope constructs a geometrical model of the scene observed along its path.The visual information is augmented by a geometrical representation,allowing the computation of geometrical parameters featuring the internal lumen of the vessel and its lesion.In case of bifurcation,multiple single meshes are merged to construct the final mesh of the whole vascular structure.A virtually active navigation system for the vascular surface reconstruction is developed on the Windows platform by Ref.[9].An example of the vascular inner surface with its centerline description is shown in Fig.2. Fig.2 Description of vascular with mesh surface and centerline1.2 ModelingFollowing the conception of multi-body representation[6],themodel is conceived as a chain of small and rigid cylindrical segments[8].Each segment represents a segment of the guide wire/catheter.It is neither compressible nor bendable and is connected to neighbors with joints,as shown in Fig.3.The black block represents the first segment which is fixed.The small spheres represent the joints and these joints can rotate through two pivots which define two degrees of freedom.The rotation of pivots is li mited by a maximum angleθmax characterizing the maximum strength.Associated with the difference in orientation between two segments is an energy measure(bending energy in the joint).The total energy of themodel combined with thevessel wall is a function of the positions of all joints.Fig.3 Multibody representation of guide wire/catheter1.3 Collision detectionThe operation of the guide wire/catheter can lead to the model collision with artery walls.A fast and interactivecollision detection algorithm is a fundamental component of the 3-D real-time simulation.Much research is addressed theissues involving the computational complexity reduction by simplifying the representation of objects in the scene[10].In this paper,the method proposed by Ref.[11]is used for the real-time collision detection,where an orthographic camera with specifications of the OpenGL graphics library is used to detect the collision.The method runs much faster than the well-known oriented boundingbox tree method[12]. The collision detection with the graphic hardware is almost independent of the artery number of polygons.Thus the real-time collision detection can be achieved even if the vasculature obtained from patient data has a great number of polygons.1.4 Cancellation of collisionIf there is a collision between the segment n and the vascular wall,the suitable position of the segment should be found to cancel the collision. Ref.[8]proposed a method based on the depth map imageanalysis,where the segment followes the direction with the maximumgradient.However,since there aremany triangular faces,a fixed length segment with a certain radius may simultaneously collide with several triangle faces.So the depth map is much morecomplex than that of the simple case,where only one plane is involved,as shown in Fig.4.It is easy tofind the direction with the maximum gradient in the simple depthmap,while it is difficult to find a suitable direction in the complex depth map.To solve the problem more efficiently,a sphere S is defined with the center point C at the end point of the segment(n-1).The radius is the segment of length(see Fig.5).All the surface points(d,O,θ)of the sphere inside thevascular can be chosen as the candidates of the end point of the segment n,the one that guarantees segment n inside the vascular and makes the angle between segments(n-1)and n minimum is chosen as the best solution.Thevasculatureis represented with many triangular faces.To decideif the segment n with the start point C is inside thevasculature,it needs to test whether the segment n intersects with any triangle face of the vasculature,since the number of triangular faces is large.the computation cost is considerably high.Fig.4 Illustration of depth mapFig.5 Collision cancellation methodTo reduce the computational cost,it firstly takes advantage of the binary CT volume data resulted from the vascular segmentation procedure.The judgment of a sphere surface point P(d,O,θ)inside the vascular becomes easy.If P is outside the vascular,the segment CP is not inside the vascular,i.e.the segment intersects with the vascular.Secondly,the bounding sphere of each triangular surface is calculated before the navigation.Only those triangular faces whose bounding spheres intersect with the sphere Sare used to test if the segment n intersects with them.Since the judgment of the sphere intersection is easy and only a few triangular faces intersect with the sphere S,so thecomputational cost is dramatically reduced.After the best point P(d,O,θ)is found,the rotational angle of th e segment n is set to be(O,θ),and the corresponding rotational angleθx andθy can be calculated.Notably,the rotation of joints is limited by the maximumangleθmax.There fore,if the rotation angleθx orθy is larger than the maximum angleθmax,an angular propag ation(AP)procedure is applied to the system[8].The procedure begins on the segment n which collides with the artery wall.Angular rotations are iteratively applied to the joint of previous segments.If thereis not suitableangle to cancel the collision,theguide wire/catheter is not suitable for the specific endovascular intervention.Then,we need to choice another one with different material properties for such a specific case.1.5 RelaxationThe above collision cancellation procedure with APcannot guarantee the whole guide model with natural status.To make it more realistic,″home springs″[13]is used to relax the guide model towards its equilibrium status after being deformed. The ″home springs″ model connects springs between the current position of each joint and the virtual zero-bending energy position[8].Thevirtual position is defined with the zero angle to the previous segments,as shown in Fig.6.Thus,the system combined with joints and segments tries to return to its equilibrium status while keep itselfinside the vascular.Please see Ref.[8]for the details of solving the constrained problem.Fig.6 ″Home springs″model2 RESULT ANALYSESTo evaluate the behavior of the proposed methodology,a system is developed on the Windows platform with VC++ and OpenGL graphic libraries,as shown is Fig.7.Given an insertion point and direction,the properties of the model by means of sliders in the left panel,the system automatically guides the motion of the multi-body model inside a specific patient vascular.A phantom(see Fig.8)and a real patient vascular(seeFig.2)segmented from CT data are tested.″Home springs″parameters used by the relaxation procedure are listed in Table 1.The evaluation is performed on a PCwith Xeon 2.33 GHz◦2 CPU and 3 GB memory.For the real patient vascular represented with 11 172 triangular meshes,the cost time for each collision detection is less than 1 ms,and that for each collision cancellation is around 50 ms, The last relaxation procedure takes less than 1 s for 50 segments.Thus,the real-time simulation can be achieved.Fig.7 Virtual reality system for 3-D multi-body simulationFig.8 Phantom with centerlineTable 1 Parameter setting of″home-springs″mi Δt k spring k collision V i1.0 0.01 15 1 500 102.1 Phantom datasetFor the phantom dataset in Fig.8,the be-haviors of the guide model withdifferent strengths and segment lengths are tested.Results before and after the relaxation are shown in Fig.9 and Table 2.The results show that the behaviors of the guide model are related to the strength e and the segment length l.Fig.9 Simulation results of guide model with different strengths and lengths of segment inside phantom2.2 Patient datasetFor the patient dataset in Fig.2,the behaviors of the guide model with different strengths and segment lengths are tested.Results before and after the relaxation are shown in Fig.10 and Table 3. The results show that different strengths and segment lengths lead to different behaviors.From Figs.9,10 and Tables 2,3,it can be seen that the guide model tends to bemore realistic and its energy becomes smaller after the relaxation.Table 2 Energy comparison before and af ter relaxation for phantom datasetParameter Energy before relaxation/J Energy after relaxation/Je=20,l=4 235.2 186.3 e=20,l=6 233.9 187.2 e=25,l=4 228.7 184.1e=25,l=6 203.5 172.3Fig.10 Simulation relaxation results of guide model with different strengths and segment lengths inside specific patient vascularTable 3 Energy comparison before and af ter relaxation for patient datasetParameter Energy before relaxation/J Energy after relaxation/Je=20,l=4 1 124.1 717.0 e=20,l=6 767.5 548.7 e=25,l=4 1 332.9 790.9e=25,l=6 779.3 526.53 CONCLUSIONThis paper presents a system for the realtime 3-D simulation of the guide wire or the catheter motion inside the specific vascular.Results show that the system is effective and promising.Theartery is segmented from CT data and represented as a 3-D mesh surface.The guide wire/catheter is modeled as a multi-body representation.The OpenGL orthographic camera is used for the real-time collision detection.While a geometry analysis combined with the AP method is developed to search the best motion direction,in this case,there is a collision.The relaxation procedure makes the simulated guide model more realistic. The behavior of the simulated guide model depends on several parameters,such as the segment length and the strength.Thus,it is necessary to do experiments to find suitable parameters for matching the physical properties of all available guide wires and catheters for the clinical usage.However,the system is still far from the goal of the endovascular intervention training.Further research is needed by considering different shapes of the guide wire tip and the catheter tip,as well as the integration of the force feedback and the active navigation algorithm in the simulation system.References：[1] Konings M K,van de Kraats E B,Alderliesten T,et al.Analytical guide wire motion algorithm for simulation of endovascular interventions[J].Med Biol Eng Comput,2003,41(6)：689-700.[2] Cotin S,Delingette H,Ayache N.Real-time elastic deformations of softtissues for surgery simulation[J].IEEE Trans Vis Comput Graph,1999(5)：62-73.[3] Delingette H.Towards realistic soft tissue modeling in medical simulation[C]//Proc of the IEEE：Special Issue on Surgical Simulation.New York,USA：IEEE Press,1998：512-523.[4] Duriez C,Cotin S,Lenoir J,et al.New approaches to catheter navigationfor interventional radiology simulation[J].Computer AidedSurgery,2006(11)：300-308.[5] Bhat S,Kesavadas T,Hoffmann K R.Aphysicallybased model for guidewire simulation on patientspecific data[J].International Congress Series,2005(1281)：479-484.[6] Cotin S,Dawson S,Meglan D,et al.ICTS,an interventional cardiology training system[C]//Proceedings of Medicine Meets Virtual Reality CA,USA：IOS Press,2000：59-65.[7] Featherstone R.The calculation of robot dynamics using axticulated-body inertias[J]. International Journal of Robotics Research,1983(2)：13-30.[8] Guilloux V,Haigron P.Simulation of guidewire navigation in complex vascular structures[C]//Proc of SPIE on Medical Imaging.Orlando：SPIE Press,2006(6141)：1-11.[9] Zhou Zhengdong,Haigron P,Shu Huazhong,etc.Optimization of intravascular brachytherapy treatment planning in peripheralarteries[J].Comput Med Imaging Graph,2007(31)：401-407.[10]Fares C,Hamam Y.Collision detection for rigid bodies：a state of theart review[C]//15th Int Conf Computer Graphics andApplications,GraphiCon′2005.Novosibirsk Akademgorodok,Russia：[s.n.],2005.[11]Lombardo JC,Cani M P,Neyret F.Real-timecollision detection for virtual surgery[C]//Proc Computer Animation′99.California,USA：IEEE Computer Society Press,1999：82-91.[12]Gottschalk S,Lin M,Manocha D.Obb-tree：ahierarchical structure for rapid inter ference detection.[C]//Proceedings of Siggraph′96.Berlin：Springer,1996：171-180.[13]LeDuc M,Payandeh S,Dill J.Toward modeling of a suturingtask[C]//Graphics Interface Conference.New York,USA：ACM Press,2003：273-279.。

Comparison of Multiobjective Evolutionary Algorithms:Empirical ResultsEckart ZitzlerDepartment of Electrical Engineering Swiss Federal Institute of T echnology 8092Zurich,Switzerlandzitzler@tik.ee.ethz.ch Kalyanmoy DebDepartment of Mechanical Engineering Indian Institute of T echnology Kanpur Kanpur,PIN208016,Indiadeb@iitk.ac.inLothar ThieleDepartment of Electrical EngineeringSwiss Federal Institute of T echnology8092Zurich,Switzerlandthiele@tik.ee.ethz.chAbstractIn this paper,we provide a systematic comparison of various evolutionary approaches tomultiobjective optimization using six carefully chosen test functions.Each test functioninvolves a particular feature that is known to cause difficulty in the evolutionary optimiza-tion process,mainly in converging to the Pareto-optimal front(e.g.,multimodality anddeception).By investigating these different problem features separately,it is possible topredict the kind of problems to which a certain technique is or is not well suited.However,in contrast to what was suspected beforehand,the experimental results indicate a hierarchyof the algorithms under consideration.Furthermore,the emerging effects are evidencethat the suggested test functions provide sufficient complexity to compare multiobjectiveoptimizers.Finally,elitism is shown to be an important factor for improving evolutionarymultiobjective search.KeywordsEvolutionary algorithms,multiobjective optimization,Pareto optimality,test functions,elitism.1MotivationEvolutionary algorithms(EAs)have become established as the method at hand for exploring the Pareto-optimal front in multiobjective optimization problems that are too complex to be solved by exact methods,such as linear programming and gradient search.This is not only because there are few alternatives for searching intractably large spaces for multiple Pareto-optimal solutions.Due to their inherent parallelism and their capability to exploit similarities of solutions by recombination,they are able to approximate the Pareto-optimal front in a single optimization run.The numerous applications and the rapidly growing interest in the area of multiobjective EAs take this fact into account.After thefirst pioneering studies on evolutionary multiobjective optimization appeared in the mid-eighties(Schaffer,1984,1985;Fourman,1985)several different EA implementa-tions were proposed in the years1991–1994(Kursawe,1991;Hajela and Lin,1992;Fonseca c2000by the Massachusetts Institute of T echnology Evolutionary Computation8(2):173-195E.Zitzler,K.Deb,and L.Thieleand Fleming,1993;Horn et al.,1994;Srinivas and Deb,1994).Later,these approaches (and variations of them)were successfully applied to various multiobjective optimization problems(Ishibuchi and Murata,1996;Cunha et al.,1997;Valenzuela-Rend´on and Uresti-Charre,1997;Fonseca and Fleming,1998;Parks and Miller,1998).In recent years,some researchers have investigated particular topics of evolutionary multiobjective search,such as convergence to the Pareto-optimal front(Van Veldhuizen and Lamont,1998a;Rudolph, 1998),niching(Obayashi et al.,1998),and elitism(Parks and Miller,1998;Obayashi et al., 1998),while others have concentrated on developing new evolutionary techniques(Lau-manns et al.,1998;Zitzler and Thiele,1999).For a thorough discussion of evolutionary algorithms for multiobjective optimization,the interested reader is referred to Fonseca and Fleming(1995),Horn(1997),Van Veldhuizen and Lamont(1998b),and Coello(1999).In spite of this variety,there is a lack of studies that compare the performance and different aspects of these approaches.Consequently,the question arises:which imple-mentations are suited to which sort of problem,and what are the specific advantages and drawbacks of different techniques?First steps in this direction have been made in both theory and practice.On the theoretical side,Fonseca and Fleming(1995)discussed the influence of differentfitness assignment strategies on the selection process.On the practical side,Zitzler and Thiele (1998,1999)used a NP-hard0/1knapsack problem to compare several multiobjective EAs. In this paper,we provide a systematic comparison of six multiobjective EAs,including a random search strategy as well as a single-objective EA using objective aggregation.The basis of this empirical study is formed by a set of well-defined,domain-independent test functions that allow the investigation of independent problem features.We thereby draw upon results presented in Deb(1999),where problem features that may make convergence of EAs to the Pareto-optimal front difficult are identified and,furthermore,methods of constructing appropriate test functions are suggested.The functions considered here cover the range of convexity,nonconvexity,discrete Pareto fronts,multimodality,deception,and biased search spaces.Hence,we are able to systematically compare the approaches based on different kinds of difficulty and to determine more exactly where certain techniques are advantageous or have trouble.In this context,we also examine further factors such as population size and elitism.The paper is structured as follows:Section2introduces key concepts of multiobjective optimization and defines the terminology used in this paper mathematically.We then give a brief overview of the multiobjective EAs under consideration with special emphasis on the differences between them.The test functions,their construction,and their choice are the subject of Section4,which is followed by a discussion about performance metrics to assess the quality of trade-off fronts.Afterwards,we present the experimental results in Section6and investigate further aspects like elitism(Section7)and population size (Section8)separately.A discussion of the results as well as future perspectives are given in Section9.2DefinitionsOptimization problems involving multiple,conflicting objectives are often approached by aggregating the objectives into a scalar function and solving the resulting single-objective optimization problem.In contrast,in this study,we are concerned withfinding a set of optimal trade-offs,the so-called Pareto-optimal set.In the following,we formalize this 174Evolutionary Computation Volume8,Number2Comparison of Multiobjective EAs well-known concept and also define the difference between local and global Pareto-optimalsets.A multiobjective search space is partially ordered in the sense that two arbitrary so-lutions are related to each other in two possible ways:either one dominates the other or neither dominates.D EFINITION1:Let us consider,without loss of generality,a multiobjective minimization problem with decision variables(parameters)and objectives:Minimizewhere(1) and where is called decision vector,parameter space,objective vector,and objective space.A decision vector is said to dominate a decision vector(also written as) if and only if(2)Additionally,in this study,we say covers()if and only if or.Based on the above relation,we can define nondominated and Pareto-optimal solutions: D EFINITION2:Let be an arbitrary decision vector.1.The decision vector is said to be nondominated regarding a set if and only if thereis no vector in which dominates;formally(3)If it is clear within the context which set is meant,we simply leave it out.2.The decision vector is Pareto-optimal if and only if is nondominated regarding.Pareto-optimal decision vectors cannot be improved in any objective without causing a degradation in at least one other objective;they represent,in our terminology,globally optimal solutions.However,analogous to single-objective optimization problems,there may also be local optima which constitute a nondominated set within a certain neighbor-hood.This corresponds to the concepts of global and local Pareto-optimal sets introduced by Deb(1999):D EFINITION3:Consider a set of decision vectors.1.The set is denoted as a local Pareto-optimal set if and only if(4)where is a corresponding distance metric and,.A slightly modified definition of local Pareto optimality is given here.Evolutionary Computation Volume8,Number2175E.Zitzler,K.Deb,and L.Thiele2.The set is called a global Pareto-optimal set if and only if(5) Note that a global Pareto-optimal set does not necessarily contain all Pareto-optimal solu-tions.If we refer to the entirety of the Pareto-optimal solutions,we simply write“Pareto-optimal set”;the corresponding set of objective vectors is denoted as“Pareto-optimal front”.3Evolutionary Multiobjective OptimizationT wo major problems must be addressed when an evolutionary algorithm is applied to multiobjective optimization:1.How to accomplishfitness assignment and selection,respectively,in order to guide thesearch towards the Pareto-optimal set.2.How to maintain a diverse population in order to prevent premature convergence andachieve a well distributed trade-off front.Often,different approaches are classified with regard to thefirst issue,where one can distinguish between criterion selection,aggregation selection,and Pareto selection(Horn, 1997).Methods performing criterion selection switch between the objectives during the selection phase.Each time an individual is chosen for reproduction,potentially a different objective will decide which member of the population will be copied into the mating pool. Aggregation selection is based on the traditional approaches to multiobjective optimization where the multiple objectives are combined into a parameterized single objective function. The parameters of the resulting function are systematically varied during the same run in order tofind a set of Pareto-optimal solutions.Finally,Pareto selection makes direct use of the dominance relation from Definition1;Goldberg(1989)was thefirst to suggest a Pareto-basedfitness assignment strategy.In this study,six of the most salient multiobjective EAs are considered,where for each of the above categories,at least one representative was chosen.Nevertheless,there are many other methods that may be considered for the comparison(cf.Van Veldhuizen and Lamont(1998b)and Coello(1999)for an overview of different evolutionary techniques): Among the class of criterion selection approaches,the Vector Evaluated Genetic Al-gorithm(VEGA)(Schaffer,1984,1985)has been chosen.Although some serious drawbacks are known(Schaffer,1985;Fonseca and Fleming,1995;Horn,1997),this algorithm has been a strong point of reference up to now.Therefore,it has been included in this investigation.The EA proposed by Hajela and Lin(1992)is based on aggregation selection in combination withfitness sharing(Goldberg and Richardson,1987),where an individual is assessed by summing up the weighted objective values.As weighted-sum aggregation appears still to be widespread due to its simplicity,Hajela and Lin’s technique has been selected to represent this class of multiobjective EAs.Pareto-based techniques seem to be most popular in thefield of evolutionary mul-tiobjective optimization(Van Veldhuizen and Lamont,1998b).In particular,the 176Evolutionary Computation Volume8,Number2Comparison of Multiobjective EAs algorithm presented by Fonseca and Fleming(1993),the Niched Pareto Genetic Algo-rithm(NPGA)(Horn and Nafpliotis,1993;Horn et al.,1994),and the Nondominated Sorting Genetic Algorithm(NSGA)(Srinivas and Deb,1994)appear to have achieved the most attention in the EA literature and have been used in various studies.Thus, they are also considered here.Furthermore,a recent elitist Pareto-based strategy,the Strength Pareto Evolutionary Algorithm(SPGA)(Zitzler and Thiele,1999),which outperformed four other multiobjective EAs on an extended0/1knapsack problem,is included in the comparison.4Test Functions for Multiobjective OptimizersDeb(1999)has identified several features that may cause difficulties for multiobjective EAs in1)converging to the Pareto-optimal front and2)maintaining diversity within the population.Concerning thefirst issue,multimodality,deception,and isolated optima are well-known problem areas in single-objective evolutionary optimization.The second issue is important in order to achieve a well distributed nondominated front.However,certain characteristics of the Pareto-optimal front may prevent an EA fromfinding diverse Pareto-optimal solutions:convexity or nonconvexity,discreteness,and nonuniformity.For each of the six problem features mentioned,a corresponding test function is constructed following the guidelines in Deb(1999).We thereby restrict ourselves to only two objectives in order to investigate the simplest casefirst.In our opinion,two objectives are sufficient to reflect essential aspects of multiobjective optimization.Moreover,we do not consider maximization or mixed minimization/maximization problems.Each of the test functions defined below is structured in the same manner and consists itself of three functions(Deb,1999,216):Minimizesubject to(6)whereThe function is a function of thefirst decision variable only,is a function of the remaining variables,and the parameters of are the function values of and.The test functions differ in these three functions as well as in the number of variables and in the values the variables may take.D EFINITION4:We introduce six test functions that follow the scheme given in Equa-tion6:The test function has a convex Pareto-optimal front:(7)where,and.The Pareto-optimal front is formed with.The test function is the nonconvex counterpart to:(8) Evolutionary Computation Volume8,Number2177E.Zitzler,K.Deb,and L.Thielewhere,and.The Pareto-optimal front is formed with.The test function represents the discreteness feature;its Pareto-optimal front consists of several noncontiguous convex parts:(9)where,and.The Pareto-optimal front is formed with.The introduction of the sine function in causes discontinuity in the Pareto-optimal front.However, there is no discontinuity in the parameter space.The test function contains local Pareto-optimal fronts and,therefore,tests for the EA’s ability to deal with multimodality:(10)where,,and.The global Pareto-optimal front is formed with,the best local Pareto-optimal front with.Note that not all local Pareto-optimal sets are distinguishable in the objective space.The test function describes a deceptive problem and distinguishes itself from the other test functions in that represents a binary string:(11)where gives the number of ones in the bit vector(unitation),ififand,,and.The true Pareto-optimal front is formed with,while the best deceptive Pareto-optimal front is represented by the solutions for which.The global Pareto-optimal front as well as the local ones are convex.The test function includes two difficulties caused by the nonuniformity of the search space:first,the Pareto-optimal solutions are nonuniformly distributed along the global Pareto front (the front is biased for solutions for which is near one);second,the density of the solutions is lowest near the Pareto-optimal front and highest away from the front:(12)where,.The Pareto-optimal front is formed with and is nonconvex.We will discuss each function in more detail in Section6,where the corresponding Pareto-optimal fronts are visualized as well(Figures1–6).178Evolutionary Computation Volume8,Number2Comparison of Multiobjective EAs5Metrics of PerformanceComparing different optimization techniques experimentally always involves the notion of performance.In the case of multiobjective optimization,the definition of quality is substantially more complex than for single-objective optimization problems,because the optimization goal itself consists of multiple objectives:The distance of the resulting nondominated set to the Pareto-optimal front should be minimized.A good(in most cases uniform)distribution of the solutions found is desirable.Theassessment of this criterion might be based on a certain distance metric.The extent of the obtained nondominated front should be maximized,i.e.,for each objective,a wide range of values should be covered by the nondominated solutions.In the literature,some attempts can be found to formalize the above definition(or parts of it)by means of quantitative metrics.Performance assessment by means of weighted-sum aggregation was introduced by Esbensen and Kuh(1996).Thereby,a set of decision vectors is evaluated regarding a given linear combination by determining the minimum weighted-sum of all corresponding objective vectors of.Based on this concept,a sample of linear combinations is chosen at random(with respect to a certain probability distribution),and the minimum weighted-sums for all linear combinations are summed up and averaged.The resulting value is taken as a measure of quality.A drawback of this metric is that only the“worst”solution determines the quality value per linear combination. Although several weight combinations are used,nonconvex regions of the trade-off surface contribute to the quality more than convex parts and may,as a consequence,dominate the performance assessment.Finally,the distribution,as well as the extent of the nondominated front,is not considered.Another interesting means of performance assessment was proposed by Fonseca and Fleming(1996).Given a set of nondominated solutions,a boundary function divides the objective space into two regions:the objective vectors for which the corre-sponding solutions are not covered by and the objective vectors for which the associated solutions are covered by.They call this particular function,which can also be seen as the locus of the family of tightest goal vectors known to be attainable,the attainment surface. T aking multiple optimization runs into account,a method is described to compute a median attainment surface by using auxiliary straight lines and sampling their intersections with the attainment surfaces obtained.As a result,the samples represented by the median attain-ment surface can be relatively assessed by means of statistical tests and,therefore,allow comparison of the performance of two or more multiobjective optimizers.A drawback of this approach is that it remains unclear how the quality difference can be expressed,i.e.,how much better one algorithm is than another.However,Fonseca and Fleming describe ways of meaningful statistical interpretation in contrast to the other studies considered here,and furthermore,their methodology seems to be well suited to visualization of the outcomes of several runs.In the context of investigations on convergence to the Pareto-optimal front,some authors(Rudolph,1998;Van Veldhuizen and Lamont,1998a)have considered the distance of a given set to the Pareto-optimal set in the same way as the function defined below.The distribution was not taken into account,because the focus was not on this Evolutionary Computation Volume8,Number2179E.Zitzler,K.Deb,and L.Thielematter.However,in comparative studies,distance alone is not sufficient for performance evaluation,since extremely differently distributed fronts may have the same distance to the Pareto-optimal front.T wo complementary metrics of performance were presented in Zitzler and Thiele (1998,1999).On one hand,the size of the dominated area in the objective space is taken under consideration;on the other hand,a pair of nondominated sets is compared by calculating the fraction of each set that is covered by the other set.The area combines all three criteria(distance,distribution,and extent)into one,and therefore,sets differing in more than one criterion may not be distinguished.The second metric is in some way similar to the comparison methodology proposed in Fonseca and Fleming(1996).It can be used to show that the outcomes of an algorithm dominate the outcomes of another algorithm, although,it does not tell how much better it is.We give its definition here,because it is used in the remainder of this paper.D EFINITION5:Let be two sets of decision vectors.The function maps the ordered pair to the interval:(13)The value means that all solutions in are dominated by or equal to solutions in.The opposite,,represents the situation when none of the solutions in are covered by the set.Note that both and have to be considered,since is not necessarily equal to.In summary,it may be said that performance metrics are hard to define and it probably will not be possible to define a single metric that allows for all criteria in a meaningful way.Along with that problem,the statistical interpretation associated with a performance comparison is rather difficult and still needs to be answered,since multiple significance tests are involved,and thus,tools from analysis of variance may be required.In this study,we have chosen a visual presentation of the results together with the application of the metric from Definition5.The reason for this is that we would like to in-vestigate1)whether test functions can adequately test specific aspects of each multiobjective algorithm and2)whether any visual hierarchy of the chosen algorithms exists.However, for a deeper investigation of some of the algorithms(which is the subject of future work), we suggest the following metrics that allow assessment of each of the criteria listed at the beginning of this section separately.D EFINITION6:Given a set of pairwise nondominating decision vectors,a neighborhood parameter(to be chosen appropriately),and a distance metric.We introduce three functions to assess the quality of regarding the parameter space:1.The function gives the average distance to the Pareto-optimal set:min(14)Recently,an alternative metric has been proposed in Zitzler(1999)in order to overcome this problem. 180Evolutionary Computation Volume8,Number2Comparison of Multiobjective EAs 2.The function takes the distribution in combination with the number of nondominatedsolutions found into account:(15) 3.The function considers the extent of the front described by:max(16) Analogously,we define three metrics,,and on the objective space.Letbe the sets of objective vectors that correspond to and,respectively,and and as before:min(17)(18)max(19)While and are intuitive,and(respectively and)need further explanation.The distribution metrics give a value within the interval()that reflects the number of-niches(-niches)in().Obviously,the higher the value,the better the distribution for an appropriate neighborhood parameter(e.g.,means that for each objective vector there is no other objective vector within-distance to it).The functions and use the maximum extent in each dimension to estimate the range to which the front spreads out.In the case of two objectives,this equals the distance of the two outer solutions.6Comparison of Different Evolutionary Approaches6.1MethodologyWe compare eight algorithms on the six proposed test functions:1.A random search algorithm.2.Fonseca and Fleming’s multiobjective EA.3.The Niched Pareto Genetic Algorithm.4.Hajela and Lin’s weighted-sum based approach.5.The Vector Evaluated Genetic Algorithm.6.The Nondominated Sorting Genetic Algorithm.Evolutionary Computation Volume8,Number2181E.Zitzler,K.Deb,and L.Thiele7.A single-objective evolutionary algorithm using weighted-sum aggregation.8.The Strength Pareto Evolutionary Algorithm.The multiobjective EAs,as well as,were executed times on each test problem, where the population was monitored for nondominated solutions,and the resulting non-dominated set was taken as the outcome of one optimization run.Here,serves as an additional point of reference and randomly generates a certain number of individuals per generation according to the rate of crossover and mutation(but neither crossover and mutation nor selection are performed).Hence,the number offitness evaluations was the same as for the EAs.In contrast,simulation runs were considered in the case of, each run optimizing towards another randomly chosen linear combination of the objec-tives.The nondominated solutions among all solutions generated in the runs form the trade-off front achieved by on a particular test function.Independent of the algorithm and the test function,each simulation run was carried out using the following parameters:Number of generations:250Population size:100Crossover rate:0.8Mutation rate:0.01Niching parameter share:0.48862Domination pressure dom:10The niching parameter was calculated using the guidelines given in Deb and Goldberg (1989)assuming the formation of ten independent niches.Since uses genotypic fitness sharing on,a different value,share,was chosen for this particular case. Concerning,the recommended value for dom of the population size wastaken(Horn and Nafpliotis,1993).Furthermore,for reasons of fairness,ran with a population size of where the external nondominated set was restricted to.Regarding the implementations of the algorithms,one chromosome was used to en-code the parameters of the corresponding test problem.Each parameter is represented by bits;the parameters only comprise bits for the deceptive function. Moreover,all approaches except were realized using binary tournament selection with replacement in order to avoid effects caused by different selection schemes.Further-more,sincefitness sharing may produce chaotic behavior in combination with tournament selection,a slightly modified method is incorporated here,named continuously updated shar-ing(Oei et al.,1991).As requires a generational selection mechanism,stochastic universal sampling was used in the implementation.6.2Simulation ResultsIn Figures1–6,the nondominated fronts achieved by the different algorithms are visualized. Per algorithm and test function,the outcomes of thefirstfive runs were unified,and then the dominated solutions were removed from the union set;the remaining points are plotted in thefigures.Also shown are the Pareto-optimal fronts(lower curves),as well as additional reference curves(upper curves).The latter curves allow a more precise evaluation of the obtained trade-off fronts and were calculated by adding max minto the values of the Pareto-optimal points.The space between Pareto-optimal and 182Evolutionary Computation Volume8,Number2f101234f2RANDFFGA NPGA HLGA VEGA NSGA SOEA SPEAFigure 1:T est function(convex).f101234f2RANDFFGA NPGA HLGA VEGA NSGA SOEA SPEAFigure 2:T est function(nonconvex).183f11234f2RANDFFGA NPGA HLGA VEGA NSGA SOEA SPEA Figure 3:T est function(discrete).f1010203040f2RANDFFGA NPGA HLGA VEGA NSGA SOEA SPEAFigure 4:T est function(multimodal).184f10246f2RANDFFGA NPGAHLGA VEGA NSGASOEA SPEA Figure 5:T est function (deceptive).f12468f2RANDFFGA NPGA HLGA VEGA NSGASOEA SPEA Figure 6:T est function(nonuniform).185reference fronts represents about of the corresponding objective space.However,the curve resulting from the deceptive function is not appropriate for our purposes,since it lies above the fronts produced by the random search algorithm.Instead,we consider all solutions with,i.e.,for which the parameters are set to the deceptive attractors (for).In addition to the graphical presentation,the different algorithms were assessed in pairs using the metric from Definition5.For an ordered algorithm pair,there is a sample of values according to the runs performed.Each value is computed on the basis of the nondominated sets achieved by and with the same initial population. Here,box plots are used to visualize the distribution of these samples(Figure7).A box plot consists of a box summarizing of the data.The upper and lower ends of the box are the upper and lower quartiles,while a thick line within the box encodes the median. Dashed appendages summarize the spread and shape of the distribution.Furthermore,the shortcut in Figure7stands for“reference set”and represents,for each test function,a set of equidistant points that are uniformly distributed on the corresponding reference curve.Generally,the simulation results prove that all multiobjective EAs do better than the random search algorithm.However,the box plots reveal that,,anddo not always cover the randomly created trade-off front completely.Furthermore,it can be observed that clearly outperforms the other nonelitist multiobjective EAs regarding both distance to the Pareto-optimal front and distribution of the nondominated solutions.This confirms the results presented in Zitzler and Thiele(1998).Furthermore, it is remarkable that performs well compared to and,although some serious drawbacks of this approach are known(Fonseca and Fleming,1995).The reason for this might be that we consider the off-line performance here in contrast to other studies that examine the on-line performance(Horn and Nafpliotis,1993;Srinivas and Deb,1994). On-line performance means that only the nondominated solutions in thefinal population are considered as the outcome,while off-line performance takes the solutions nondominated among all solutions generated during the entire optimization run into account.Finally,the best performance is provided by,which makes explicit use of the concept of elitism. Apart from,it even outperforms in spite of substantially lower computational effort and although uses an elitist strategy as well.This observation leads to the question of whether elitism would increase the performance of the other multiobjective EAs.We will investigate this matter in the next section.Considering the different problem features separately,convexity seems to cause the least amount of difficulty for the multiobjective EAs.All algorithms evolved reasonably distributed fronts,although there was a difference in the distance to the Pareto-optimal set.On the nonconvex test function,however,,,and have difficulties finding intermediate solutions,as linear combinations of the objectives tend to prefer solutions strong in at least one objective(Fonseca and Fleming,1995,4).Pareto-based algorithms have advantages here,but only and evolved a sufficient number of nondominated solutions.In the case of(discreteness),and are superior to both and.While the fronts achieved by the former cover about of the reference set on average,the latter come up with coverage.Among the considered test functions,and seem to be the hardest problems,since none of the algorithms was able to evolve a global Pareto-optimal set.The results on the multimodal problem indicateNote that outside values are not plotted in Figure7in order to prevent overloading of the presentation. 186。

Matlab的第三方工具箱大全(按住CTRL点击连接就可以到达每个工具箱的主页面来下载了)Matlab Toolboxes∙ADCPtools - acoustic doppler current profiler data processing∙AFDesign - designing analog and digital filters∙AIRES - automatic integration of reusable embedded software∙Air-Sea - air-sea flux estimates in oceanography∙Animation - developing scientific animations∙ARfit - estimation of parameters and eigenmodes of multivariate autoregressive methods∙ARMASA - power spectrum estimation∙AR-Toolkit - computer vision tracking∙Auditory - auditory models∙b4m - interval arithmetic∙Bayes Net - inference and learning for directed graphical models∙Binaural Modeling - calculating binaural cross-correlograms of sound∙Bode Step - design of control systems with maximized feedback∙Bootstrap - for resampling, hypothesis testing and confidence interval estimation ∙BrainStorm - MEG and EEG data visualization and processing∙BSTEX - equation viewer∙CALFEM - interactive program for teaching the finite element method∙Calibr - for calibrating CCD cameras∙Camera Calibration∙Captain - non-stationary time series analysis and forecasting∙CHMMBOX - for coupled hidden Markov modeling using max imum likelihood EM ∙Classification - supervised and unsupervised classification algorithms∙CLOSID∙Cluster - for analysis of Gaussian mixture models for data set clustering∙Clustering - cluster analysis∙ClusterPack - cluster analysis∙COLEA - speech analysis∙CompEcon - solving problems in economics and finance∙Complex - for estimating temporal and spatial signal complexities∙Computational Statistics∙Coral - seismic waveform analysis∙DACE - kriging approximations to computer models∙DAIHM - data assimilation in hydrological and hydrodynamic models∙Data Visualization∙DBT - radar array processing∙DDE-BIFTOOL - bifurcation analysis of delay differential equations∙Denoise - for removing noise from signals∙DiffMan - solv ing differential equations on manifolds∙Dimensional Analysis -∙DIPimage - scientific image processing∙Direct - Laplace transform inversion via the direct integration method∙DirectSD - analysis and design of computer controlled systems with process-oriented models∙DMsuite - differentiation matrix suite∙DMTTEQ - design and test time domain equalizer design methods∙DrawFilt - drawing digital and analog filters∙DSFWAV - spline interpolation with Dean wave solutions∙DWT - discrete wavelet transforms∙EasyKrig∙Econometrics∙EEGLAB∙EigTool - graphical tool for nonsymmetric eigenproblems∙EMSC - separating light scattering and absorbance by extended multiplicative signal correction∙Engineering Vibration∙FastICA - fixed-point algorithm for ICA and projection pursuit∙FDC - flight dynamics and control∙FDtools - fractional delay filter design∙FlexICA - for independent components analysis∙FMBPC - fuzzy model-based predictive control∙ForWaRD - Fourier-wavelet regularized deconvolution∙FracLab - fractal analysis for signal processing∙FSBOX - stepwise forward and backward selection of features using linear regression∙GABLE - geometric algebra tutorial∙GAOT - genetic algorithm optimization∙Garch - estimating and diagnosing heteroskedasticity in time series models∙GCE Data - managing, analyzing and displaying data and metadata stored using the GCE data structure specification∙GCSV - growing cell structure visualization∙GEMANOVA - fitting multilinear ANOVA models∙Genetic Algorithm∙Geodetic - geodetic calculations∙GHSOM - growing hierarchical self-organizing map∙glmlab - general linear models∙GPIB - wrapper for GPIB library from National Instrument∙GTM - generative topographic mapping, a model for density modeling and data visualization∙GVF - gradient vector flow for finding 3-D object boundaries∙HFRadarmap - converts HF radar data from radial current vectors to total vectors ∙HFRC - importing, processing and manipulating HF radar data∙Hilbert - Hilbert transform by the rational eigenfunction expansion method∙HMM - hidden Markov models∙HMMBOX - for hidden Markov modeling using maximum likelihood EM∙HUTear - auditory modeling∙ICALAB - signal and image processing using ICA and higher order statistics∙Imputation - analysis of incomplete datasets∙IPEM - perception based musical analysisJMatLink - Matlab Java classesKalman - Bayesian Kalman filterKalman Filter - filtering, smoothing and parameter estimation (using EM) for linear dynamical systemsKALMTOOL - state estimation of nonlinear systemsKautz - Kautz filter designKrigingLDestimate - estimation of scaling exponentsLDPC - low density parity check codesLISQ - wavelet lifting scheme on quincunx gridsLKER - Laguerre kernel estimation toolLMAM-OLMAM - Levenberg Marquardt with Adaptive Momentum algorithm for training feedforward neural networksLow-Field NMR - for exponential fitting, phase correction of quadrature data and slicing LPSVM - Newton method for LP support vector machine for machine learning problems LSDPTOOL - robust control system design using the loop shaping design procedure LS-SVMlabLSVM - Lagrangian support vector machine for machine learning problemsLyngby - functional neuroimagingMARBOX - for multivariate autogressive modeling and cross-spectral estimation MatArray - analysis of microarray dataMatrix Computation- constructing test matrices, computing matrix factorizations, visualizing matrices, and direct search optimizationMCAT - Monte Carlo analysisMDP - Markov decision processesMESHPART - graph and mesh partioning methodsMILES - maximum likelihood fitting using ordinary least squares algorithmsMIMO - multidimensional code synthesisMissing - functions for handling missing data valuesM_Map - geographic mapping toolsMODCONS - multi-objective control system designMOEA - multi-objective evolutionary algorithmsMS - estimation of multiscaling exponentsMultiblock - analysis and regression on several data blocks simultaneously Multiscale Shape AnalysisMusic Analysis - feature extraction from raw audio signals for content-based music retrievalMWM - multifractal wavelet modelNetCDFNetlab - neural network algorithmsNiDAQ - data acquisition using the NiDAQ libraryNEDM - nonlinear economic dynamic modelsNMM - numerical methods in Matlab textNNCTRL - design and simulation of control systems based on neural networks NNSYSID - neural net based identification of nonlinear dynamic systemsNSVM - newton support vector machine for solv ing machine learning problems NURBS - non-uniform rational B-splinesN-way - analysis of multiway data with multilinear modelsOpenFEM - finite element developmentPCNN - pulse coupled neural networksPeruna - signal processing and analysisPhiVis- probabilistic hierarchical interactive visualization, i.e. functions for visual analysis of multivariate continuous dataPlanar Manipulator - simulation of n-DOF planar manipulatorsPRT ools - pattern recognitionpsignifit - testing hyptheses about psychometric functionsPSVM - proximal support vector machine for solving machine learning problems Psychophysics - vision researchPyrTools - multi-scale image processingRBF - radial basis function neural networksRBN - simulation of synchronous and asynchronous random boolean networks ReBEL - sigma-point Kalman filtersRegression - basic multivariate data analysis and regressionRegularization ToolsRegularization Tools XPRestore ToolsRobot - robotics functions, e.g. kinematics, dynamics and trajectory generation Robust Calibration - robust calibration in statsRRMT - rainfall-runoff modellingSAM - structure and motionSchwarz-Christoffel - computation of conformal maps to polygonally bounded regions SDH - smoothed data histogramSeaGrid - orthogonal grid makerSEA-MAT - oceanographic analysisSLS - sparse least squaresSolvOpt - solver for local optimization problemsSOM - self-organizing mapSOSTOOLS - solving sums of squares (SOS) optimization problemsSpatial and Geometric AnalysisSpatial RegressionSpatial StatisticsSpectral MethodsSPM - statistical parametric mappingSSVM - smooth support vector machine for solving machine learning problems STATBAG - for linear regression, feature selection, generation of data, and significance testingStatBox - statistical routinesStatistical Pattern Recognition - pattern recognition methodsStixbox - statisticsSVM - implements support vector machinesSVM ClassifierSymbolic Robot DynamicsTEMPLAR - wavelet-based template learning and pattern classificationTextClust - model-based document clusteringTextureSynth - analyzing and synthesizing visual texturesTfMin - continous 3-D minimum time orbit transfer around EarthTime-Frequency - analyzing non-stationary signals using time-frequency distributions Tree-Ring - tasks in tree-ring analysisTSA - uni- and multivariate, stationary and non-stationary time series analysisTSTOOL - nonlinear time series analysisT_Tide - harmonic analysis of tidesUTVtools - computing and modifying rank-revealing URV and UTV decompositions Uvi_Wave - wavelet analysisvarimax - orthogonal rotation of EOFsVBHMM - variation Bayesian hidden Markov modelsVBMFA - variational Bayesian mixtures of factor analyzersVMT- VRML Molecule Toolbox, for animating results from molecular dynamics experimentsVOICEBOXVRMLplot - generates interactive VRML 2.0 graphs and animationsVSVtools - computing and modifying symmetric rank-revealing decompositions WAFO - wave analysis for fatique and oceanographyWarpTB - frequency-warped signal processingWAVEKIT - wavelet analysisWaveLab - wavelet analysisWeeks - Laplace transform inversion via the Weeks methodWetCDF - NetCDF interfaceWHMT - wavelet-domain hidden Markov tree modelsWInHD - Wavelet-based inverse halftoning via deconvolutionWSCT - weighted sequences clustering toolkitXMLTree - XML parserYAADA - analyze single particle mass spectrum dataZMAP - quantitative seismicity analysis。

用MATLAB的时重要MATLAB toolsADCPtools - acoustic doppler current profiler data processingAFDesign - designing analog and digital filtersAIRES - automatic integration of reusable embedded softwareAir-Sea - air-sea flux estimates in oceanographyAnimation - developing scientific animationsARfit - estimation of parameters and eigenmodes of multivariate autoregressive methodsARMASA - power spectrum estimationAR-Toolkit - computer vision trackingAuditory - auditory modelsb4m - interval arithmeticBayes Net - inference and learning for directed graphical modelsBinaural Modeling - calculating binaural cross-correlograms of soundBode Step - design of control systems with maximized feedbackBootstrap - for resampling, hypothesis testing and confidence interval estimationBrainStorm - MEG and EEG data visualization and processingBSTEX - equation viewerCALFEM - interactive program for teaching the finite element methodCalibr - for calibrating CCD camerasCamera CalibrationCaptain - non-stationary time series analysis and forecastingCHMMBOX - for coupled hidden Markov modeling using maximum likelihood EM Classification - supervised and unsupervised classification algorithmsCLOSIDCluster - for analysis of Gaussian mixture models for data set clusteringClustering - cluster analysisClusterPack - cluster analysisCOLEA - speech analysisCompEcon - solving problems in economics and financeComplex - for estimating temporal and spatial signal complexitiesComputational StatisticsCoral - seismic waveform analysisDACE - kriging approximations to computer modelsDAIHM - data assimilation in hydrological and hydrodynamic modelsData VisualizationDBT - radar array processingDDE-BIFTOOL - bifurcation analysis of delay differential equationsDenoise - for removing noise from signalsDiffMan - solving differential equations on manifoldsDimensional Analysis -DIPimage - scientific image processingDirect - Laplace transform inversion via the direct integration methodDirectSD - analysis and design of computer controlled systems with process-oriented modelsDMsuite - differentiation matrix suiteDMTTEQ - design and test time domain equalizer design methodsDrawFilt - drawing digital and analog filtersDSFWAV - spline interpolation with Dean wave solutionsDWT - discrete wavelet transformsEasyKrigEconometricsEEGLABEigTool - graphical tool for nonsymmetric eigenproblemsEMSC - separating light scattering and absorbance by extended multiplicative signal correctionEngineering VibrationFastICA - fixed-point algorithm for ICA and projection pursuitFDC - flight dynamics and controlFDtools - fractional delay filter designFlexICA - for independent components analysisFMBPC - fuzzy model-based predictive controlForWaRD - Fourier-wavelet regularized deconvolutionFracLab - fractal analysis for signal processingFSBOX - stepwise forward and backward selection of features using linear regressionGABLE - geometric algebra tutorialGAOT - genetic algorithm optimizationGarch - estimating and diagnosing heteroskedasticity in time series modelsGCE Data - managing, analyzing and displaying data and metadata stored using the GCE data structure specificationGCSV - growing cell structure visualizationGEMANOVA - fitting multilinear ANOVA modelsGenetic AlgorithmGeodetic - geodetic calculationsGHSOM - growing hierarchical self-organizing mapglmlab - general linear modelsGPIB - wrapper for GPIB library from National InstrumentGTM - generative topographic mapping, a model for density modeling and data visualizationGVF - gradient vector flow for finding 3-D object boundariesHFRadarmap - converts HF radar data from radial current vectors to total vectorsHFRC - importing, processing and manipulating HF radar dataHilbert - Hilbert transform by the rational eigenfunction expansion methodHMM - hidden Markov modelsHMMBOX - for hidden Markov modeling using maximum likelihood EMHUTear - auditory modelingICALAB - signal and image processing using ICA and higher order statisticsImputation - analysis of incomplete datasetsIPEM - perception based musical analysisJMatLink - Matlab Java classesKalman - Bayesian Kalman filterKalman Filter - filtering, smoothing and parameter estimation (using EM) for linear dynamical systemsKALMTOOL - state estimation of nonlinear systemsKautz - Kautz filter designKrigingLDestimate - estimation of scaling exponentsLDPC - low density parity check codesLISQ - wavelet lifting scheme on quincunx gridsLKER - Laguerre kernel estimation toolLMAM-OLMAM - Levenberg Marquardt with Adaptive Momentum algorithm for training feedforward neural networksLow-Field NMR - for exponential fitting, phase correction of quadrature data and slicingLPSVM - Newton method for LP support vector machine for machine learning problemsLSDPTOOL - robust control system design using the loop shaping design procedure LS-SVMlabLSVM - Lagrangian support vector machine for machine learning problemsLyngby - functional neuroimagingMARBOX - for multivariate autogressive modeling and cross-spectral estimationMatArray - analysis of microarray dataMatrix Computation - constructing test matrices, computing matrix factorizations, visualizing matrices, and direct search optimizationMCAT - Monte Carlo analysisMDP - Markov decision processesMESHPART - graph and mesh partioning methodsMILES - maximum likelihood fitting using ordinary least squares algorithmsMIMO - multidimensional code synthesisMissing - functions for handling missing data valuesM_Map - geographic mapping toolsMODCONS - multi-objective control system designMOEA - multi-objective evolutionary algorithmsMS - estimation of multiscaling exponentsMultiblock - analysis and regression on several data blocks simultaneouslyMultiscale Shape AnalysisMusic Analysis - feature extraction from raw audio signals for content-based music retrievalMWM - multifractal wavelet modelNetCDFNetlab - neural network algorithmsNiDAQ - data acquisition using the NiDAQ libraryNEDM - nonlinear economic dynamic modelsNMM - numerical methods in Matlab textNNCTRL - design and simulation of control systems based on neural networksNNSYSID - neural net based identification of nonlinear dynamic systemsNSVM - newton support vector machine for solving machine learning problemsNURBS - non-uniform rational B-splinesN-way - analysis of multiway data with multilinear modelsOpenFEM - finite element developmentPCNN - pulse coupled neural networksPeruna - signal processing and analysisPhiVis - probabilistic hierarchical interactive visualization, i.e. functions for visual analysis of multivariate continuous dataPlanar Manipulator - simulation of n-DOF planar manipulatorsPRTools - pattern recognitionpsignifit - testing hyptheses about psychometric functionsPSVM - proximal support vector machine for solving machine learning problemsPsychophysics - vision researchPyrTools - multi-scale image processingRBF - radial basis function neural networksRBN - simulation of synchronous and asynchronous random boolean networksReBEL - sigma-point Kalman filtersRegression - basic multivariate data analysis and regressionRegularization ToolsRegularization Tools XPRestore ToolsRobot - robotics functions, e.g. kinematics, dynamics and trajectory generationRobust Calibration - robust calibration in statsRRMT - rainfall-runoff modellingSAM - structure and motionSchwarz-Christoffel - computation of conformal maps to polygonally bounded regions SDH - smoothed data histogramSeaGrid - orthogonal grid makerSEA-MAT - oceanographic analysisSLS - sparse least squaresSolvOpt - solver for local optimization problemsSOM - self-organizing mapSOSTOOLS - solving sums of squares (SOS) optimization problemsSpatial and Geometric AnalysisSpatial RegressionSpatial StatisticsSpectral MethodsSPM - statistical parametric mappingSSVM - smooth support vector machine for solving machine learning problemsSTATBAG - for linear regression, feature selection, generation of data, and significance testingStatBox - statistical routinesStatistical Pattern Recognition - pattern recognition methodsStixbox - statisticsSVM - implements support vector machinesSVM ClassifierSymbolic Robot DynamicsTEMPLAR - wavelet-based template learning and pattern classificationTextClust - model-based document clusteringTextureSynth - analyzing and synthesizing visual texturesTfMin - continous 3-D minimum time orbit transfer around EarthTime-Frequency - analyzing non-stationary signals using time-frequency distributions Tree-Ring - tasks in tree-ring analysisTSA - uni- and multivariate, stationary and non-stationary time series analysisTSTOOL - nonlinear time series analysisT_Tide - harmonic analysis of tidesUTVtools - computing and modifying rank-revealing URV and UTV decompositions Uvi_Wave - wavelet analysisvarimax - orthogonal rotation of EOFsVBHMM - variation Bayesian hidden Markov modelsVBMFA - variational Bayesian mixtures of factor analyzersVMT - VRML Molecule Toolbox, for animating results from molecular dynamics experimentsVOICEBOXVRMLplot - generates interactive VRML 2.0 graphs and animationsVSVtools - computing and modifying symmetric rank-revealing decompositionsWAFO - wave analysis for fatique and oceanographyWarpTB - frequency-warped signal processingWAVEKIT - wavelet analysisWaveLab - wavelet analysisWeeks - Laplace transform inversion via the Weeks methodWetCDF - NetCDF interfaceWHMT - wavelet-domain hidden Markov tree modelsWInHD - Wavelet-based inverse halftoning via deconvolutionWSCT - weighted sequences clustering toolkitXMLTree - XML parserYAADA - analyze single particle mass spectrum dataZMAP - quantitative seismicity analysis。

Frobenius MethodFrobenius method is a powerful mathematical technique used to find solutions to differential equations, especially those with singular points. It was developed by Ferdinand Georg Frobenius, a German mathematician, in the late 19th century. The method has applications in various fields of science and engineering, including physics, chemistry, and biology.IntroductionDifferential equations are mathematical equations that involve derivatives of unknown functions. They are used to model a wide range of phenomena in the physical and natural sciences. Solving these equations is crucial in understanding the behavior of systems governed by differential equations.The Frobenius method is particularly useful when dealing with second-order linear differential equations that have regular singular points. A regular singular point is a point at which the coefficients of the differential equation become infinite or indeterminate.ProcedureTo apply the Frobenius method, we follow these steps:1.Assume a power series solution: We assume that the solution to thedifferential equation can be expressed as an infinite series ofpowers of x.2.Substitute the power series into the differential equation: Wesubstitute the assumed power series into the original differential equation and equate coefficients of like powers of x.3.Determine recurrence relations: By equating coefficients, weobtain recurrence relations that relate each coefficient toprevious coefficients in the power series.4.Solve for initial values: Using initial conditions or boundaryconditions, we can solve for some or all of the initial values in terms of known quantities.5.Determine convergence: We analyze the radius and interval ofconvergence for the power series solution to ensure its validitywithin a certain domain.6.Construct the general solution: Once we have determined allnecessary coefficients, we can construct the general solution as a linear combination of terms involving these coefficients.ExampleLet’s consider an example to illustrate how to apply the Frob enius method:We want to solve Bessel’s equation:x^2 y’’ + x y’ + (x^2 - n^2) y = 0where n is a constant.1.Assume a power series solution:We assume that the solution can be expressed as a power series:y(x) = Σ(a_n * x^(n+r))where r is an arbitrary constant.2.Substitute the power series into the differential equation: Substituting the power series into Bessel’s equation, we obtain:Σ[(n+r)(n+r-1) * a_n * x^(n+r-2) + (n+r) * a_n * x^(n+r-1) + (x^2 - n^2) * a_n * x^(n+r)] = 03.Determine recurrence relations:By equating coefficients of like powers of x, we obtain the following recurrence relation:a_(n+2) = -((n+r)^2 - n^2)/(n+1)(n+2) * a_n4.Solve for initial values:To solve for initial values, we need to consider the behavior of the coefficients as n approaches negative infinity. We can chooseappropriate values for r to ensure convergence and determine all necessary initial values.5.Determine convergence:We analyze the recurrence relation to determine the radius and interval of convergence for the power series solution.6.Construct the general solution:Once we have determined all necessary coefficients, we can construct the general solution by combining terms involving these coefficients.ApplicationsThe Frobenius method has various applications in science and engineering. Some examples include:1.Quantum mechanics: The Frobenius method is used to solveSchrödinger’s equation in atomic physics and quantum mechanics.2.Heat conduction: The Frobenius method can be applied to solvepartial differential equations governing heat conduction invarious materials.3.Fluid dynamics: The Frobenius method is used to solve differentialequations describing fluid flow and turbulence in fluid dynamics.4.Electromagnetism: The Frobenius method is used to solve Maxwell’sequations, which describe the behavior of electric and magneticfields.ConclusionThe Frobenius method provides a systematic approach for findingsolutions to differential equations with regular singular points. By assuming a power series solution and equating coefficients, we can determine the necessary recurrence relations and initial values to construct the general solution. This method has wide-rangingapplications in various scientific and engineering disciplines, makingit an essential tool for researchers and practitioners alike.。