Lyapunov指数（精选3篇）
Lyapunov指数 篇1
实时准确的空中交通复杂性评估是解决未来空中交通高密度运行的有效工具。近年来, 随着空中交通流量的不断增长, 空中交通复杂性研究日益受到关注。当前, 空中交通复杂性研究的主要方法是NASA提出的动态密度概念, 它是基于交通几何特征、冲突几何特征、协调动作及交通结构等多重变量在内的线性加权函数。由于动态密度权重确定带有极强的主观因素, 且其应用的动力学模型经过极大简化, 未考虑非线性因素, 因此不能满足未来空中交通复杂性评估的需求[1]。针对动态密度存在的问题, 国外陆续有人提出基于内禀属性的空中交通复杂性研究[2,3], 虽然避免了动态密度评价复杂度时的主观权重, 但由于采用构建解析模型的方法, 对数据样本的选择过于依赖, 且解析模型过于线性化, 依然不能较好地反应空中交通的复杂本质。
空中交通系统, 作为一个复杂的非线性系统, 其复杂程度可用系统可预测性的大小来表征, 可预测性强表示系统复杂程度低, 反之可预测性低表示系统复杂程度高, 即定量表征为系统的Lyapunov指数, 以可预测程度作为空中交通系统复杂性指标, 研究空中交通复杂度, 有助于理解空中交通复杂度产生的原因和其内在演化趋势。
针对扇区内飞机对时序距离, 利用时序分析的方法, 根据重构时序数据与空中交通系统数学模型几何等价这一关系, 判断空中交通系统可预测程度。采用Lyapunov指数为指标, 识别飞机间的态势, 建立系统可预测程度、飞机态势、复杂度三者间的联系, 利用厦门扇区雷达数据验证Lyapunov指数识别空中交通态势的有效性。
本文首先介绍复杂动力系统时序分析的流程, 然后给出Lyapunov指数计算方法, 最后分别给出正常状态和事故征候下的应用实例。
1 复杂动力系统时序分析
在空中交通系统中飞机可以被视为运动的粒子, 由粒子组成的动力系统中, 向量场决定了动力系统的复杂性。由于空中交通系统高度复杂, 通常只能通过观测手段获得外部单变量时间序列。而时间序列分析是对系统的观测数据进行分析, 找出空中交通系统的特征、内在变化规律从而揭示空中交通系统复杂度的过程。
对于未知数学模型的空中交通系统, 可以通过雷达得到其对应的向量场
针对向量场中位置子空间的时间序列{xn}进行研究, 做以下相空间重构。
从而形成m维状态空间, 在重构的m维状态空间中可以建立数学模型
由Takens定理可知, 只要选取合适的m和τ, 就可以利用时序数据重构与原未知数学模型的动力系统的几何特征等价的相空间。对于空中交通系统的解释、分析都可以在重构的相空间下进行[4,5]。过程如图1所示。
2 空中交通态势识别算法
Lyapunov指数是表征相空间内轨迹分离速率的量, 用于度量非线性时间序列的可预测程度。系统的Lyapunov指数表征随时间变化时, 系统对于初值的敏感程度, Lyapunov指数为负, 说明该系统在对应的相空间中体积在该方向上是收缩的, 而在该方向上的运动是稳定的;相反Lyapunov指数为正, 说明系统在对应相空间中的体积在该方向上是膨胀的, 也就说明在该方向上的运动越来越不相关, 成为长期不可预测的运动。
对于空中交通而言, 短期稳定可预测的运动状态就可以对应较简单态势, 而长期发散不可预测的运动状态就可以对应复杂态势。态势的简单与复杂, 可直接对应复杂度。
2.1 相空间重构
重构相空间是计算Lyapunov指数的前提, 相空间重构主要是为确定平均时间间隔τ与嵌入维数m。利用自相关函数法求时间延时, 利用预测误差最小法计算嵌入维, 从而将空中交通系统进行相空间重构。
2.1.1 自相关函数求时间延迟
由于雷达数据长度有限, 且存在噪声污染与测量误差, 因此时间延迟τ不能任意选择。现通过自相关函数Cl (τ) , 判断雷达数据xn和xn+τ是否线性相关, 从而求得时间延迟τ。
雷达数据的自相关函数如下。
当Cl (τ) 首次为零时, 在平均意义下雷达数据xn与xn+τ线性无关, 此时τ为所求时间延迟[6]。
2.1.2 预测误差最小法计算嵌入维
对于雷达数据重构的相空间xn, 利用临近轨道欧氏距离的最大值范数之比的平均值E (m) , 计算嵌入维数m。
由于E (m) 的变化只与嵌入维m和延迟时间τ间隔相关, 为求得m增大过程中的最小m, 重新定义E1 (m) 为
当m大于某一个m0时, E1 (m) 趋于平稳, 则m0+1即为所求最小嵌入维[7]。
2.2 计算最大Lyapunov指数
根据上节所述的算法, 对雷达数据{xn}Nn=1进行重构。
以x1为基准点, 在相点集{xn}, n=1, 2, …, N- (m-1) τ中可找到与x1欧式距离最近的一点, 并求得两点间距离L1。经过适当的时间步长T1演化后, L1变为L'1, 可得指数增长率为
在演化过后的基准点附近, 取与演化后基准点距离最小点, 求得距离L2, 继续经过时间步长T2演化后得到指数增长率为
遍历空间内所有点, 记演化的总步数为M, 则取指数增长率λi的作为雷达数据最大Lyapunov指数的估计值[8], 即
3 算例
3.1 扇区交通流时间序列
为识别扇区空中交通态势, 在厦门扇区流量高峰时间段某日9:30~9:45, 进行数据采集, 采样以时间和厦门一号扇区为共同限制条件, 采样数据为各飞机雷达航迹数据。将飞机雷达航迹数据转化为各飞机对之间的距离, 由此就得到了实测飞机对距离时序数据。这些将作为识别扇区空中交通态势的原始依据。
3.2 飞机对间混沌判别
利用第2节所述算法, 计算各飞机对间时序距离的Lyapunov指数。结果如下。
由计算结果可知, 飞机对2、3、4、7、10, 其时序距离的Lyapunov指数为负, 不存在可预测性, 而飞机对1、4、5、8、9, 其时序距离的Lyapunov指数为正, 存在可预测性。
为更好地体现态势与Lyapunov指数之间的关系, 针对某一事故征候的雷达数据, 将雷达数据转化为时序距离数据, 计算时序距离数据的Lyapunov指数。
计算出两机之间时序距离的最大Lyapunov指数为0.268 8。此时飞机对间距离的Lyapunov指数值较大, 说明系统不可预测的程度较大, 此事故征候的态势为两机侧向间隔超标。
3.3 空中交通复杂态势的识别与产生的原因
观察计算结果可以发现, 最大Lyapunov指数分为两类:正、负, 分别对应跟驰、反向、大角度汇聚、对头基本四种态势。
同一航线上同航向飞行的飞机流, 呈现跟驰态势, 前后机之间的距离主要是纵向间隔, 保持正常巡航速度且前后机之间纵向间隔较大, 此时对飞机下一时刻的态势预判较为容易, 所以lyapunov指数为负, 且飞机对2、6、7、10对应的Lyapunov指数绝对值较大, 说明跟驰态势的有序性最强。
反向态势从飞机对进入扇区开始两机距离不断增大, 两机间互不影响, 因此也不需要对下一时刻进行预判, 但由于开始距离较近, 且距离变化趋势上可能存在波动, 飞机对3的Lyapunov指数虽然小于零, 但绝对值极小, 可能由于某些因素的波动在将来时间段, 转化为不可预测状态, 使空中交通复杂度增加。
随着扇区流量的增大, 会出现各个航向的飞机流, 此时扇区内飞机对间, 可能存在大角度汇聚、对头飞行等态势, 管制员为避免冲突, 不断调整飞机速度、航向、高度, 这个过程中飞机对间距离不断变化, 从而导致飞机对间的不可预测状态。观察飞机对Lyapunov值大小, 可以发现大角度汇聚态势飞机对1、4的值较大而对头态势飞机对5、7、9值较小, 说明大角度汇聚的不可预测程度明显高于对头飞行的不可预测程度。
随着流量的增加, 飞机流的有序运动将变为无序的不可预测状态。在无序的不可预测状态中, 受到外界某些因素的干扰, 可能改变系统的可预测性。研究表明, 飞机流间存在可预测性的不同, 空中交通态势可以用可预测程度进行识别。管制员能力差异、飞机性能差异和空域结构等原因, 导致了混沌现象的产生。在实际管制工作中, 认为对头飞行与大角度会聚的复杂性要高于同向跟驰的情况。而间隔小于标准时, 认为此时复杂度大大提高, 导致管制员不能很好地对飞机之间间隔进行调配。Lyapunov指数态势识别算法的结果与实际管制工作可较好吻合, 该算法可不依赖于具体数学模型, 利用客观雷达数据对飞机间态势进行识别。
4 结束语
给出了利用雷达数据计算飞机对复合间隔Lyapunov指数的过程与公式。首次讨论了基于飞机对距离的Lyapunov指数识别飞机间态势, 在空中交通中的应用, 并通过实际数据算例和实际管制情况的对照, 证明了该方法的可行性和有效性。
由于扇区结构、飞机速度等原因导致雷达数据有限, 但并不影响研究结果的一般性与普遍性, 随着空中交通监视手段的不断进步, 可以获得更短时更高频的数据, 以期分析影响Lyapunov指数的具体因素, 从而实现精细化的空中交通态势识别与空中交通复杂度控制。
参考文献
[1] Laudeman I, Shelden S, Branstrom R.Dynamic density:an air traffic mangement metric.National Technical Information Service, 1998
[2] Delahaye D, Puechmorel S.Air traffic complexity:towards intrinsic metrics.3rd USA/Europe ATM R&D Seminar, 2000
[3] Delahaye D, Paimblanc P, Puechmorel S, et al.A new air traffic complexity metric based on dynamical system modelization.21st Digital Avionics Systems Conference, 2002;1:4A2-1—4A2-12
[4] 王海燕, 卢山.非线性时间序列分析及其应用.北京:科学出版社, 2002:12—21Wang Haiyan, Lu Shan.Nonlinear time series analysis and its application.Beijing:Science Press, 2002:12—21
[5] Takens F.Detecting strange attractors in turbulence.Berlin:Lecture Notes in Mathematics, 1981:366—381
[6] Abarbanel H D I, Brown R, Sidorowich J J, et al.The analysis of observed data in physical system.Reviews of Modern Physics, 1993;65 (4) :1331—1392
[7] Cao L Y.Practical method for determining the minimum embedding dimension of a scalar time series.Physica D, 1997;110 (5) :43 —50
[8] Wolf A.Determing Lyapunov exponents from a time series.Physica D, 1985;16 (3) :285—317
Lyapunov指数 篇2
变形监测一个重要目的在于防患未然,即检查各种工程建筑物和地质构造的稳定性,对变形进行分析,掌握规律,准确进行预测预报,从而为有关部门决策提供依据。由于变形体变形机理的复杂性和多样性,对变形分析与建模理论和方法的研究,需引入系统科学及非线性科学的理论,采用数学模型来逼近、模拟和揭示变形体的变形规律和动态特征。关于预测预报的方法很多,但都有一定的局限性,没有一种方法能够完全适应各种情况下的预测和分析,因此需要不断地改进和完善。
混沌作为非线性科学的一个重要分支,是当今倍受科学界关注的前沿学科和研究热点。20世纪90年代以来,混沌科学与其他科学相互渗透,在很多领域得到了广泛的应用。在混沌的应用上,根据混沌系统提取的非线形时间序列对系统的未来进行预测是研究的重点。
本文将混沌时间序列理论引入对沉降监测数据的预测分析,建立了基于Lyapunov指数的变形预测模型,通过工程实例,对该模型进行了验证,得到一些有益的结论,能够对混沌时间序列在变形监测中更广泛的应用提供参考价值。
1 Lyapunov指数
混沌是一种貌似无规则的运动,是在确定性非线性系统中不需要附加任何随机因素出现类似随机行为。混沌系统的最大特点在于系统的演化对初始条件十分敏感。
混沌的离散情况常常表现为混沌时间序列,混沌时间序列是由混沌模型生成的具有混沌特性的时间序列,混沌时间序列中蕴涵着系统丰富的动力学信息,混沌时间序列是混沌理论通向现实世界的一个桥梁,是混沌的一个重要应用领域。Lyapunov指数是量化初始闭轨道的指数发散和估计系统混沌量,它从整体上反映了动力系统的混沌量水平,因此,基于混沌时间序列的Lyapunov指数计算和预测显得尤其重要。
在一维动力系统an+1=f(an)中,初始两点迭代后是互相分离的,还是靠拢的,关键取决于的值。若,则迭代使得两点分开,若,则迭代使得两点靠拢。但是在不断的迭代过程中,的值也分随之发生变化,分离与靠拢交替进行。为了表示从整体上看相邻两状态分离的情况,必须对时间(迭代次数)取平均。因此,设平均每次迭代所引起的指数中的指数为λ于是原来相距ψ的两点经过n次迭代的相距为
取极限ψ→0,n→∞,式(1)变为
式(2)通过变形计算可简化为
式(3)中的λ与初始值的选取有很大关系,称为原动力系统的李雅普诺夫(Lyapunov)指数,它表示系统在多次迭代中平均每次迭代所引起的指数分离中的指数。
对一个耗散动力系统而言,判断其是否为混沌的重要标志就是看该系统的最大Lyapunov指数是否为正。Lyapunov指数为正说明了相邻运动轨迹在相空间中的发散趋势,是“初值敏感性”(SIC)特性的定量化体现。SIC意味着初始时刻一个微小的波动会导致系统将来极大的变化,也就是说,系统相空间中初始相邻的轨迹会伴随着系统演化而迅速分离,分离的速率可以通过最大Lyapunov指数λ刻画。由于Lyapunov指数λ刻画的是系统全局范围的发散性质,同时其具体数值不因对原系统的线性变换而改变,因此可作为系统混沌特性的可靠度量。在已知系统观测序列Sn的情况下,λ的大小可直接从序列Sn中估计出来。
2 基于最大Lyapunov指数变形预测模型的建立
2.1 重构相空间
相空间是一个以系统变量数为维数的多维空间。相空间中的一个点(即相点),代表了动力系统在某一时刻的一个特定状态。相空间里相点的连线,构成了点在相空间的轨道,即相轨道。相轨道表示了系统状态随时间的演变。
根据相空间重构理论,x(t 1),x(t 2),...,x(t n)是所研究的时间序列,可得m维延迟矢量:
其中τ称为延迟时间,m称为嵌入维数。建立了一个多维相空间Y i(i=1,2,...,n),得到一组时间序列x(ti),相空间中的点的个数n=N-(m-1)τ。
2.2 利用C-C方法计算延迟时间和嵌入维数
C-C方法是于1999年由D.Kugiumtzis,R.Eykholt和J.D.Salas提出,该方法应用关联积分能够同时估计出τd和τω。时间延迟τd确保xi个成分相互依赖,但不依赖于m;而时间窗口τω依赖于m,且τ随m而变化
C-C方法计算步骤如下:
(1)计算给定的沉降时间序列的标准差σ。
(2)编程计算下列三个量:
分别如下:
C(m,r,t)为关联积分的极限值[1]。
(3)依据上述计算结果:
的第一个极小值t对应时间延迟τd=tτs;
的第一个极小值t对应时间延迟τd=tτs;
,最小值t对应时间窗口τω=tτs。
τs为变形监测时间序列的监测周期。
2.3 Wolf方法计算最大Lyapunov指数
从单变量的时间序列提取Lyapunov指数的方法仍然是基于时间序列的重构相空间。Wolf等人(1985)提出直接基于相轨平面、相体积等演化来估计Lyapunov指数。这类方法统称为Wolf方法,它在混沌的基础研究和基于Lyapunov指数的混沌时间序列预测中应用十分广泛。本文就是采用Wolf方法计算Lyapunov指数。
设混沌时间序列x1,x2,...,xk,...嵌入维数m,时间延迟τ,则重构相空间
取初始点Y(t 0)设其与最近邻点Y0(t 0)的距离为0L,追踪这两点的时间演化,直到1t时刻,其间距超过某规定值,保留Y(t 1),并在Y(1t)邻近另找一个点Y1(t 1),使得,并且与之夹角尽可能的小,继续上述过程,直至Y(t)到达时间序列的终点N,这时追踪演化过程总的迭代次数为M,则最大Lyapunov指数为
2.4基于Lyapunov指数的预测模式
Lyapunov指数可用来表征沉降变形系统的混沌行为,沉降变形系统在相空间中相邻轨迹的指数发散,刻画相空间中相体积收缩和膨胀的几何特性,是一个很好的预报沉降结果的参数。
不妨设YM为预报的中心点,相空间中YM的最近的邻点为Yk,其距离为d M(0),最大Lyapunov指数为λ1,即
其中点YM+1只有最后一个分量x(tn+1)未知,故x(tn+1)是可预报的。式(12)就是基于Lyapunov指数的预测模式,图1是基于Lyapunov指数的变形预测流程图。
3 工程应用实例
为检验基于Lyapunov指数的变形预测模型的有效性,采用樊村隧道的沉降监测数据进行验证。该隧道位于浙江省杭(州)-金(华)-衢(州)高速公路衢州窑上段第B合同段,隧道起止桩号K253+280~+525,其中明洞30m(进口端20m,出口端10m),暗洞215m,隧道几何线形和净空按100km/h设计,隧道照明部分按80km/h设计。另外,该隧道位于侵蚀丘陵区,地势起伏较大,沿隧道纵向分布的主要岩体有:含碎石亚粘土、含粘性士碎石、泥质粉砂岩、粉砂岩、石英质砂岩、含砂砾岩、砂砾石等不同岩性的岩体,地质情况极为复杂。为了避免开挖隧道而出现危害,在2002年4月25日起,利用隧道位移实时监测系统对该隧道拱顶下沉进行监测,一天三次,每次约半小时,本次采用的典型点D5的监测数据,此点在监测期间共采集了355期沉降监测数据。
首先确定时间延迟和嵌入维,根据C-C方法原理,使用Matlab编程计算求出时间延迟τd和嵌入维m。,结果如图2所示。
然后使用Wolf方法计算Lyapunov指数λ1,并对第326~355天的沉降值进行预测,预测结果如图3所示。
为了检验预测结果的好坏,我们定义相对误差为:
相对误差=|预测值-真实值|/真实值
得出预测结果的相对误差值图,如图4所示。从图3和图4可以看出,基于Lyapunov指数的变形预测模型在此工程实例中,预测结果跟实际监测结果基本相符,相对误差在0~0.09之间波动,故判断该预测模型的预测结果精度很好。
4结论
经过对监测的沉降时间序列进行相空间重构,使用基于最大Lyapunov指数的预测方法来预测隧道拱顶的沉降,得到较高精度的预测结果,同时可得如下结论:
(1)计算的Lyapunov指数表明,隧道拱顶的沉降具有混沌性质。
(2)Lyapunov指数预测结果表明,通过选择合适的嵌入维数,能够增加预测的准确程度。
(3)由于混沌的初值的敏感性,使用Lyapunov指数进行长期预测不可能,但只要有足够好的模型和对初始条件的精确观察,它的确定性在预测能力消失之前可以进行短期预测。
摘要：变形危害巨大,变形监测与变形预测则成为必然,由于变形的过程受到地质、水文、地震和人类工程活动等因素的影响,可视为一种具有混沌特征的动力系统。故本文以混沌理论为基础,提出了基于最大Lyapunov指数的变形观测模式,并利用工程实例进行了分析研究。研究表明,变形监测时间序列中的最大Lyapunov指数均大于0,根据混沌理论,可判断变形序列存在混沌现象,同时预测的结果显示,基于混沌时间序列的最大Lyapunov指数法的预测具有较高的精度,故研究成果具有的一定的理论和应用价值,为变形预测提供了新途径。
关键词：混沌,时间序列,变形,Lyapunov指数,预测
参考文献
[1]吕金虎,陆君安,陈士华.混沌时间序列分析及其应用[M].武汉:武汉大学出版社,2002.
[2]蒋廷臣,张勤等.基于小波方法的非线性回归模型研究[J].测绘学报,2006,35(4):337~341.
[3]蒋传文,候志俭,张勇传.基于关联度的高嵌入维混沌预测方法研究[J].系统工程与电子技术,2002,24(12):65~66,103.
[4]鲁铁定,周世健,张立亭,官云兰.测量平差教学中MATLAB软件的应用[J].地矿测绘,2004,20(1):43~45.
[5]张步涵,刘小华,万建平等.基于混沌时间序列的负荷预测及其关键问题的分析[J].电网技术,2004,28(13):32~49.
[6]梁勇,孟桥,陆佶人.Lyapunov指数的算法改进与加权预测[J].声学技术,2006,25(5):463~469.
[7]CHEN Baohua1,LI Jianping,DING Ruiqiang.Nonlinear local Lyapunov exponent and atmospheric predictability research.Science in China Series D:Earth Sciences,2006,49(10)1111~1120.
[8]董昭.带跳的非线性随机微分方程的Lyapunov指数的估计[J].应用数数学学报,2007,30(1):23~36.
[9]冯明库,丘水生,晋建秀.一种Lyapunov指数算法及其实现[J].计算机应用,2007,27(1):50~54.
Lyapunov指数 篇3
对于方程已知的动力学系统,可以很方便的计算出整个系统的Lyapunov指数谱[2]。对于时间序列的Lyapunov指数,到目前为止,计算的方法很多,但它们大体上分属于两类:Wolf方法[3]和Jocobian方法[4]。但这两种方法实际操作起来都比较困难,Wolf方法适用于时间序列无噪声,切空间中小向量的演变高度非线性;Jocobian方法适用于时间序列噪声大,切空间中小向量的演变接近线性。1993年Rsenstein et al.[5]提出了计算时间序列最大Lyapunov指数的小数据量方法,该方法操作起来比较方便,而且具有以下优点:(1)该方法对小数据量可靠;(2)计算量并不大;(3)相对容易操作。但当数据量增大时,该方法的计算量很大。本文将Delaunay三角剖分[5]应用到该算法中,通过理论分析和对几种离散映射的仿真实验表明,该算法大大减小了计算的时间,而且对嵌入维数不敏感。
1 基本原理
1.1 时间序列最大Lyapunov指数的小数据量算法
在实际中,测得的往往只是单一变量的时间序列,而系统任一分量的演化是由与之相互作用着的其它分量所决定的,所以单一变量的时间序列隐含着整个系统的运动规律。因此,首先进行相空间重构。设混沌时间序列{x1,x2,…,xN},选取适当的嵌入维数m和时间延迟τ,则重构相空间:
在重构相空间后,找出相空间中每个点Yj的最近邻点Yj^,然后对相空间中的每个点Yj,计算出该邻点对在第i个离散时间步后的距离dj(i),对每个i,求出所有j的lndj(i)平均y(i),即
(2)式中q是非零dj(i)的数目,并用最小二乘法做出回归直线,该直线的斜率就是最大的Lyapunov指数[5]。
1.2 Delaunay三角剖分在Lyapunov指数计算中的应用
平面上给定n个点p1,p2,…,pn,所谓平面点集三角剖分是指用互不相交直线段连接pi与pj;1≤i,j≤n,i≠j,并使凸壳内的每个区域是一个三角形。1934年俄国科学家Delaunay[6]提出了一种特殊的三角剖分———Delaunay三角剖分,是最近点意义下的Voronoi图的直线对偶。Delaunay三角剖分在计算几何中的应用非常广泛。其中一个最重要应用就是寻找最近邻点。
在1.1所描述的计算时间序列最大Lyapunov指数的算法中,在相空间重构后,需要查找每个点的最近邻点。因此,可以将Delaunay三角剖分应用到Lyapunov指数的计算中。具体步骤如下:对于相空间重构后的每个点构成的点集,首先,利用Guibas-Stolfi算法[7]构造出点集的Delaunay三角剖分。有定理[8]指出,在点集的delaunay三角剖分中,每个点和它的最近邻点之间都存在一条边。然后,对于点集中的每个点pi,查找与它之间有边存在的所有点,找出其最近邻点。
1.3 算法时间复杂度分析
在整个Lyapunov指数的计算过程中,主要的计算量都集中在寻找重构相空间之后每个点的最近邻点中。如果我们按照常规的对每个点依次查询的方法,那么,假设数据的长度为N,找出一个点的最近邻点需查询(N-1)次,那么N个点需要查询N×(N-1)次,计算的时间复杂度为O(N[2])。由此可见,当数据长度较长时,计算的时间复杂度很大,即需要的计算时间很长。
那么用Delaunay三角剖分来寻找最近邻点,计算的复杂度又如何呢?Delaunay三角剖分是最近点意义下的Voronoi图的对偶图。因此,首先构造点集的Voronoi图,再做其直线对偶图,即对Voronoi图的每条边作通过点集中某两点的垂线,便得到Delaunay三角剖分。构造Voronoi图所需的运算次数,在定理[9]中指出:用O(Nlg N)时间能够造平面上N个点的一个集合的Vonoroi图,而且是最优的。那么,如果数据长度为N,构造Delaunay三角剖分的计算量为O(Nlg N)。然后根据Delaunay三角剖分查询最近邻点。因为任何三角剖分最多有(3N-6)条边[9],而每条边为两个点所共有,因此,最多需要查询2×(3N-6)次,查找的时间复杂度为O(N)。那么总的时间复杂度为O(Nlg N)+O(N)=O(Nlg N)。而前已指出,按照常规方法查询每个点的最近邻点所需的时间为O(N[2]),可见,特别是当数据量较大时,本文提出的算法在时间复杂度上有较大的提高。
2 仿真实验
根据前面的算法讨论,取下面两组数据对新算法进行验证,并分析实验结果。
2.1 实验1 Logistic映射
这是描述生物种群系统演化的典型模型,即虫口模型,其中xi表示亲代数量,xi+1表示子代(约化)数量,μ为常数,这里取μ=4.0。重构相空间时,选取时间延迟为1,嵌入维数为3。所用的数据点为5 000。由理论推导[10]可知,最大Lyapunov指数的理论值为λ=ln(2)=0.693 15。用该算法计算出的最大Lyapunov指数为0.687 88。算法时间复杂度的比较:用本文的方法计算时间为25.562 s,而用通常的计算邻点的方法计算时间为939.969 s。同时,改变相空间重构中的嵌入维数m,分别计算最大Lyapunov指数(如图1),并计算误差(见表1)。
注:虚线为最小二乘拟合后的直线,该直线的斜率为最大Lyapunov指数。
2.2 实验2 Hénon映射
这里设a=1.4,b=0.3。相空间重构时,选取时间延迟为1,嵌入维数为4。所用数据点为5 000。最大Lyapunov指数的理论值[3]为0.418。实验所得的值为0.429 77。算法时间复杂度的比较:用本文的方法计算时间为21.313 s,而用通常的寻找最近邻点的方法,计算时间为985.031 s。嵌入维数和最大Lyapunov指数之间的关系以及计算误差见表2。
3 结语
通过对Logistic映射和Hénon映射的仿真实验可以看出,计算所得的最大Lyapunov指数和理论结果之间存在一定的误差,但误差很小,最大不超过3%(见表1、表2)。而且无论对于Logistic映射还是Hénon映射,不同的嵌入维数m,计算出的最大Lyapunov指数有明显的区别。由误差分析可知:(a)当嵌入维数为2时,实验结果最接近真实值;(b)并不是嵌入维数越小,计算所得的最大Lyapunov指数值越准确,如在Hénon映射中嵌入维数为5时的最大Lyapunov指数比嵌入维数为4时更接近于理论值。但Lyapunov指数是刻画系统邻近轨道平均发散程度的量,应该与嵌入维数m的选取关系不大。而之所以得出上述试验结果,是因为在上述过程中,时间延迟τ和嵌入维数m的选取是相互独立的,而实际上它们之间应该有一定的联系。
本文提出了利用Delaunay三角剖分计算时间序列最大Lyapunov指数的算法,分析了算法的时间复杂度,并将算法应用于几种离散映射,分别计算了在不同嵌入维数下的最大Lyapunov指数。仿真结果表明,该方法减小了计算的时间复杂度,节约了计算时间;同时,嵌入维数的选择对计算结果有一定的影响。然而,在计算的精度上,由于需要线性拟合,曲线线性段的选取缺乏客观的标准,从而影响了计算的精确,这是该方法最大的缺陷,有待于进一步改进。
摘要：针对时间序列最大Lyapunov指数计算速度慢的缺陷,研究了小数据量算法,提出了基于Delaunay三角剖分的最大Lyapunov指数的计算方法。利用Delaunay三角剖分方法解决了邻点搜索速度慢的问题。详细地介绍了算法步骤,分析了算法的运算量,并应用于几种离散映射。仿真试验表明:该方法较稳定、可靠,同时对相空间重构中的嵌入维数不敏感。
关键词：Delaunay三角剖分,最大Lyapunov指数,相空间重构,混沌
参考文献
[1]刘秉正.非线性动力学与混沌基础(修订本).长春:东北师范大学出版社,1995
[2]Von Bremen H,Udwadia F.An efficient QR based method for the computation of Lyapunov exponents.Physica D,1997;101:1—16
[3]Wolf A,Swift J B.Determining Lyapunov exponents from a time se-ries.Physica D,1985;16:285—317
[4]Barana G,Tsuda I.A new method for computing Lyapunov expo-nents.Phys Lett A,1993;175:421—427
[5]Rosenstein M T,Collins J J,De Luca C J.A ptactical method for calculating largest Lyapunov exponents from small data sets.Physica D,1993;65:117—134
[6]Delaunay B.Sur la sph啨re vide.Bull Acad Sci USSR(Ⅶ)Classe Sci Mat Nat,1934;7(6):793—800
[7]Guibas L,Stolfi J.Primitives for the manipulation of general subdivi-sions and computation of Voronoi diagrams.ACM Transactions on Graphics,1985;4:74—123
[8]周培德.计算几何.北京:清华大学出版社,1999
[9]Preparata F R,Shamos M I.Computational geometry—an introduc-tion.New York:Springer-Verlag,1985