混沌行为的判据

混沌行为的判据（精选4篇）

混沌行为的判据 篇1

伴随着激光器的诞生,人们在实验上观察到激光器输出存在尖峰效应以及类似随机噪声现象,即激光器混沌现象。近年来,随着混沌理论研究的不断深入,人们逐渐意识到混沌可以被广泛地应用到各个领域中,而在光电领域的应用尤为突出。例如,利用混沌控制技术,消除激光器的混沌,实现对周期状态,尤其是高周期态的稳定控制,如激光加工(切割)与激光测距等。另外,以光学混沌作为信息载体,通过适当编码、解码过程进行保密通信。尤其是利用混沌对初始条件极其敏感的特性,将微弱信号加载到初始条件上,通过系统吸引子的变化,探测出被测信号的微小变化,达到利用激光混沌进行弱信号检测的目的。由此可见,研究激光混沌具有重要的实际意义。鉴于激光混沌的广泛应用,迄今,它已经成为光学领域研究的重要课题[1,2,3,4]。

能够产生混沌信号的激光器主要有NH3激光器、Nd:YAG激光器、CO2激光器、He-Ne激光器以及半导体激光器等[5,6,7,8,9]。基于激光器的动力学方程,运用稳定性理论可以判定激光系统的动力学行为。值得注意的是,在理论分析的过程中,一些合理的近似可能导致一些动力学行为的丢失。基于此,对均匀加宽激光器的动力学方程进行了修正。考虑光场随空间演化的效应,得到修正后的均匀加宽激光器的动力学方程。进一步计算系统的最大Lyapunov指数随参量的演化关系,发现存在Lyapunov指数大于零的参数区域。在此区域内,系统呈现丰富的时空混沌行为。

1 激光器动力学方程的修正

以单模均匀加宽NH3激光器来说明激光混沌的产生机理,其实验装置如图1所示。

图1采用三面反射镜构成环形腔,而衍射光栅采用它的一级耦合泵浦光。将泵浦光引入光学谐振腔,用它的零级衍射做反射镜。用N2O激光器的P(13)支谱线的激光做泵浦源。通过空间滤波器和小孔组合,改变泵浦光功率。在泵浦光的作用下,NH3分子的能级之间形成粒子数反转分布,从而产生激光。描述该激光器的动力学方程为[10]

其中,E′、P′、D′分别代表电场强度、宏观电极化强度和反转粒子数密度。并且,用参量γc代表激光器腔内光场的衰减速率;γ⊥表示宏观极化强度的迟豫速率;γΣ表示原子能级上粒子数衰减速率。D0表示激光器的激发源非相干激发产生原子能级上粒子数的稳定值。参量α则由下列关系式表示

其中,N′为工作物质的原子数密度;μ是工作物质原子跃迁电偶极矩的模数;ωc为谐振腔的共振频率;ε0为真空中的介电常数。

当式(1)左端为零时,E′、P′、D′将不随时间变化。此时,可以确定其稳态解Es、Ps、Ds为下列形式

作标度变换。

在这种变换下,并利用式(3)、式(1)可以变换为下列形式

当方程中的参量为σ=10,b=8/3,c>24.74,系统处于混沌态。

值得注意的是,方程(4)略去了光场的空间变化,即其中的状态变量仅随时间变化。实际上,光场随空间变化的效应是存在的。在这种情况下,方程(4)中的状态变量应该是时间和空间的函数。考虑到光场随空间变化的效应,在方程(4)中加入了光场的空间扩散项,将方程(4)修正为

其中,t仍表示时间;r表示空间坐标;d为空间扩散系数,∇2=∂2/∂r2。

基于稳定性理论可以判定系统的动力学行为。其中,计算系统的Lyapunov指数来判定系统的动力学行为是常用的有效方法之一。这里,数值计算了修正后的方程(5)的最大Lyapunov指数。取参量σ=10,b=8/3,空间扩散系数d=0.01,作出系统的最大Lyapunov指数随参量c的演化关系如图2所示。图2显示最大Lyapunov指数既存在大于零的区域,也存在小于零的区域。在参量c=28处所对应的最大Lyapunov指数大于零,表明在该参量下系统处于时空混沌态。此时系统状态变量的时空演化如图3~图5所示。

由图3~图5可以看出,方程中的状态变量随时间和空间的演化虽然呈现随机的行为,但它们又是有界的。并且改变初始值,重复模拟图3~图5,即状态变量随时间和空间的演化,发现演化规律完全不同,但仍保持随机有界特征,说明系统对初始条件的微小变化具有高度的敏感依赖性。激光混沌的这些特性在实际中的应用是非常广泛的,例如,利用其进行信息通信;利用多个激光器进行远程通信时中继信号的同步转换和逐级放大以及利用激光混沌及其特性进行图像加密等。因此,研究激光混沌及其特性具有重要的科学意义和广泛的应用价值。

2 结论

对均匀加宽激光器的动力学方程进行了修正。首先考虑光场不仅随时间变化,也随空间变化的因素,对均匀加宽激光器的动力学方程进行修正。进一步依据稳定性理论,通过计算系统的最大Lyapunov指数确定系统的动力学行为。研究结果表明,修正后的激光器动力学方程的最大Lyapunov指数存在大于零和小于零的区域。在最大Lyapunov指数大于零的区域所对应的参数之下,系统可以呈现时空混沌行为。

参考文献

[1]尹元昭.混沌同步和混沌保密通讯的实验研究[J].电子科学学报,1998,20(1):93-96.

[2]孙军强,黄德修,易河清.基于半导体激光放大器增益饱和的波长转换的小信号分析[J].物理学报,1997,46(12):2369-2375.

[3]户永清.超混沌LC振荡电路的混合反馈同步与保密通信[J].通信技术,2010,43(10):113-116.

[4]冯秀琴,姚治海,田作林,等.简并光学参量振荡器的超混沌控制与周期态同步[J].物理学报,2010,59(12):8414-8419.

[5]Roy R,Thornburg K S Jr.Experimental synchronization of chaotic lasers[J].Phys Rev Lett,1994,72(13):2009-2012.

[6]Uchida A,Higa K,Shiba T,et al.Generalized synchronization of chaos in He-Ne lasers[J].Phys Rev E,2003,68(1):016215-7.

[7]Vicente R,Tang S,Mulet J,et al.Synchronization properties of two self-oscillating semiconductor lasers subject to delayed optoelectronic mutual coupling[J].Phys Rev E,2006,73(4):047201-4.

[8]Rogister F,Roy R.Spatiotemporal chaos synchronization in a large array of locally coupled nonidentical lasers[J].Laser Phys,2005,15(2):313-319.

[9]Aviad Y,Reidler I,Kinze W,et al.Phase synchronization in mutually coupled chaotic diode lasers[J].Phys Rev E,2008,78(2):025204-4.

[10]Haken H.Analogy between higher instabilities in fluids and lasers[J].Phys Lett A,1975,53(1):77-78.

混沌行为的判据 篇2

关键词：Logistic映射,混沌,李亚普诺夫指数

引言

一个形式上非常简单的系统其动力学行为可能非常复杂,这也正是混沌系统的重要特征之一,本文将通过对Logistic映射的分析研究,揭示混沌产生的动力学机制,并揭示混沌现象中普遍成立的规律。

Verhulst建立了生物种群的繁衍模型,其主要思想是:假设某一类种群第n代的个体总数为Nn,生态环境能够维持的种群个体总数为N,假设。由于存在繁衍关系,第n+1代的个体总数必然与xn和Nn有关,定义种群的繁殖增长率为,可以得到:,如果再考虑到生物种群个体总数受环境资源的制约,第n+1代能够获取的资源必然与(1-xn)成正比,因此,xn+1必然正比与(1-xn),

(1)式被称为虫口模型,也称为单参数的Logistic映射模型[1]。线性项ax代表虫口数的平均增长率,而非线性项-ax2(a>0)体现环境资源对种群繁衍的制约因素。通过设定初值x0并研究数值序列x1,x2,…,xn,…的变化规律,我们就得到了种群繁衍的规律,然而,x1,x2,…,xn,…数值序列的渐进行为却并不简单。计算发现,Logistic映射的渐进行为与a的取值密切相关。

1 Logistic映射动力学行为

当a<1时,xn的不动点xn=0是稳定的,经过几次迭代逐渐趋向xn=0,意味着物种最终走向灭绝(如图1所示)。当1<a<3时,xn的不动点是稳定的,经过几次迭代逐渐趋向,物种最终保持数量的稳定(如图2所示)。

特别需要指出的是,存在一个特殊的a值即a∞=3.569945672,当a>a∞时,Logistic映射进入混沌区域,如图3所示。上述迭代过程意味着形式上非常简单的Logistic映射有着复杂的运动学行为,为了揭示数值序列x1,x2,…,xn,…的渐进行为的复杂性,我们可以从下面系列图示得到合理解释。当a>3时,xn的不动点进入失稳状态,经过几次迭代逐渐进入周期轨道,出现倍周期分叉现象[2],并最终进入混沌区域[3]。

2 Logistic映射混沌行为的动力学原因

图4是周期三窗口及阵发性混沌示意图,周期三存在的原因是:

当a接近周期三窗口时,系统慢慢进现阵发的混沌,每两次阵发混沌之间的近周期行为的时间长短是随机的[3],这一行为被称为阵发混沌。阵发混沌的临界点a=ac≈3.8284,当a接近ac时周期三轨道的持续时间越来越长,a→ac时周期时间趋于无穷长。a>ac时阵发混沌消失,重新回到周期三轨道[4],如图5(b)所示。在ac点上发生的分叉行为叫作切分叉。切分叉的名称来源于图5(c)上所示的分叉点出现xn+3=xn与xn+3=f(3)(xn)两线是相切的。

3 李亚普诺夫指数

假设两轨道的初始相邻

则

其中γ为李亚普诺夫指数。图6是Logistic映射的Lyapunov指数图。如果γ>0,则该系统具有初值敏感性,系统具有混沌特征。

李亚普诺夫指数的定义式为:

虽然混沌系统的动力学方程是确定的,混沌运动的长期演化结果却不确定[5]。对于已知的xn可以由(1)式精确计算出xn+1,而所谓精确计算是指我们把xn的测量精确到,则能保证xn+1的预言精度。但研究混沌运动的长期行为时我们并不能准确预言它的轨道,导致了运动结果的不确定性。当Lyapunov指数λ>0时,△xn的取值范围随n以平均eλ的速率指数增加,对xn预测的误差逐级扩大,xn最终取什么值就完全不可预测了。

4 结语

Logistic映射的混沌行为可以利用数学迭代的演化规律进行判定,并通过Lyapunov指数进行表征,当Lyapunov指数取正值时,其迭代行为就会表现为混沌。

参考文献

[1]H.O.Peitgen,H.Jurgens,D.Saupe Chaos and fractals New frontiers of science[M],384-387,2008

[2]徐渟.混沌在工程中的应用.物理与工程,No.5,53-55,2000

[3]邓绍江,李传东.混沌理论及其在密码学中的应用[J].重庆建筑大学学报,25(5),123-12,2003

[4]OttE,grebogiC,YorkeJA.Chaos/xaocSovit-AmericanPerspectivesinNonlinearScience[M].Ed.byD.k.Camp2bell,1989AmericanInst.ofPhys.NY,1990.0,153-172

混沌行为的判据 篇3

滚动轴承支承的转子系统由于轴承非线性刚度而出现复杂的非线性动力学特性。近年来, 刚性转子滚动轴承系统的分叉及混沌特性研究得到了广泛的重视。

Jang等[1]考虑滚动轴承的5 个自由度以及表面波纹度建立了刚性转子滚动轴承系统的非线性动力学方程, 发现套圈的波纹会导致径向位移和角位移响应峰值附近出现边频带。Harsha[2-3]考虑轴承间隙、不平衡力研究了转子系统的动力学响应, 发现系统的动力学行为与滚动体的通过频率有关, 且当通过频率及其谐波与固有频率相等、外圈波纹度阶数与滚动体数目相等时振动幅值将变得极大。Bai等[4]研究了主轴-滚动轴承系统非线性动力学行为, 发现轴承间隙的减小有利于提高主轴轴承系统的稳定性。高尚晗等[5]研究了主轴- 滚动轴承在负游隙情况下的机床主轴-滚动轴承系统的非线性动力学特性, 揭示了主轴系统的混沌演化过程。崔立等[6]研究了圆柱滚子轴承刚性转子系统周期运动分岔特性, 发现随着径向间隙、阻尼和力矩的变化, 周期运动将产生倍周期或Hopf分岔, 分岔转速随参数变化而改变。

以上的研究对象均为刚性转子, 模型涉及转轴的弯曲变形和陀螺力矩等参数。随着旋转机械转速的提高, 柔性转子系统的设计与分析变得越来越重要。近年也有一些研究柔性转子轴承系统的文章发表, 如:Laha等[7]研究了油膜轴承支承的柔性转子系统的分叉行为, 分析了转轴的材料、刚度和质量等参数对转子系统分叉行为的影响。Villa等[8]、Sinou[9]研究了球轴承支承的柔性转子系统的非线性动力学行为, 发现系统响应中存在跳跃现象和超谐波, 但其研究仅考虑了4个自由度的转轴节点和2个轴承自由度, 难以满足实际工况的需求。

本文采用12 自由度Euler-Bernoulli杆单元建立柔性转子滚子轴承系统的非线性动力学模型, 研究系统的混沌行为, 分析轴承结构参数与转轴结构参数对转子系统混沌行为的影响规律。

1 计算模型

图1所示为柔性转子轴承系统模型, 采用12自由度的Euler-Bernoulli杆单元, 圆盘简化为质点并考虑其质量与转动惯量, 转子由两个滚动轴承支承。

1.1 柔性转子轴承系统动力学模型

柔性转子系统的动力学方程为

式中, M为包含转轴、圆盘和轴承的总质量矩阵;C为总阻尼矩阵;G为总陀螺矩阵;K为总刚度矩阵;X为转子各节点的位移向量;f (X, t) 为包括转子各节点轴承力、重力、转子不平衡力和外载荷的矩阵。

先对转子系统进行轴段与节点划分, 然后进行转子系统的质量矩阵、刚度矩阵、阻尼矩阵、陀螺矩阵、载荷矩阵求解, 并按照节点的顺序对各矩阵进行组装, 具体过程如下:首先计算杆单元的质量矩阵、陀螺矩阵、刚度矩阵、载荷矩阵并组装;然后将刚性圆盘的质量矩阵、陀螺矩阵、不平衡力矩阵叠加到所在节点的相应矩阵中;之后将支承轴承的质量矩阵、非线性轴承力叠加到所在节点的相应矩阵中;最后计算系统的结构阻尼、轴承阻尼矩阵并组装。

1.2杆单元

图2所示为针对空间杆单元分析获得的空间杆单元两端节点位移。

针对图2 中12 自由度Euler-Bernoulli杆单元, 若每个节点考虑6个自由度, 则共有6个广义位移和6个广义力, 其表达式为

式中, xi、yi、zi分别为节点i沿x、y、z方向的线位移;θxi、θyi、θzi分别为节点i处截面绕x、y、z轴的转角位移;Fxi、Fyi、Fzi分别为节点i在x、y、z方向受到的轴向力和剪切力;Mxi、Myi、Mzi分别为节点i在x、y、z方向所受的扭矩和弯矩。

假设已知杆单元横截面面积、截面惯性矩、单元的扭转惯性矩、长度、材料弹性模量和剪切模量, 则可根据有限元理论[10]求出杆单元的质量矩阵Ms、刚度矩阵Ks、阻尼矩阵Cs、陀螺矩阵Gs。

1.3 刚性圆盘模型

当轴上安装有圆盘时, 将其视为刚性圆盘, 并将其质量矩阵、陀螺矩阵、不平衡力矩阵叠加到所在节点的相应矩阵中。

假设已知圆盘的质量、半径、不平衡质径积, 则可建立刚性圆盘的质量矩阵Md、陀螺矩阵Gd、不平衡力矩阵Fd, 其中Fd为

式中, me为圆盘的不平衡质径积;ω 为转速。

1.4 滚子轴承模型

考虑普遍受载的圆柱滚子轴承, 其模型如图3所示, 图中, Dr为滚子直径, D1、D2分别为外圈外径和内圈内径, δ0为径向间隙, le为带凸度滚子的长度, ls为带凸度滚子直线部分的长度。

假设滚子数目为N, 使用切片法将滚子分成nr个圆片。对第j个滚子进行受力分析, 假设第j个滚子方位角为j, 根据赫兹接触理论, 并使用拟动力学方法建立滚子轴承的非线性平衡方程组, 使用Newton-Raphson法可求出滚子与套圈的接触力[11]。将各滚子的接触力分解, 即可求出圆柱滚子轴承的非线性轴承力矩阵:

式中, F2jk为第j个滚子的第k个切片与滚子轴承内圈的作用力。

1.5 结构阻尼求解

Rayleigh提出的结构阻尼模型计算表达式为

式中, ε1、ε2为两个振型的阻尼比, 根据经验取ε1=0.005, ε2=0.01;ω1、ω2为转子系统的二阶固有频率。

式 (7) 中转子系统的固有频率ω1、ω2可根据计算得到的质量矩阵、刚度矩阵解|K - Mω2|= 0得到。

2 计算方法

采用Runge-Kutta法、Newton-Raphson法进行非线性动力学方程组求解, 根据FPA修正法确定求解周期, 然后求解最大Lyapunov指数, 判断系统的动力学行为。

2.1 求解周期的确定

转子轴承系统中存在轴承变刚度激励, 还可能存在不平衡力激励。不平衡力产生的激励周期为轴转动周期的整数倍, 但轴承的变刚度激励周期往往不是轴转动周期的整数倍, 所以在判断和求解时, 采用修正的FPA法建立统一的求解周期[12], 其表达式定义为

式中, T为求解周期;Td为不平衡力激励周期;TVC为轴承变刚度激励周期;ε为常数, 取ε= 0.01;K为比例系数, K =1, 2, …, nk。

根据式 (9) 进行循环计算, 直至找到满足其要求的K值, 代入式 (8) 即得求解周期。

2.2 Lyapunov指数计算

Lyapunov指数表示在相平面中2条相邻轨线间的距离随时间的平均指数发散率, 它明确地区分了确定性运动和混沌运动[13]。

对于连续系统有

设 (τ) 为一基准轨线, X (τ) 为其相邻轨线。定义矩阵X = [x1x2… xn]T的范数为

考察X (τ) 与 (τ) 之间的距离随时间延续的发散程度, 即

设在τ0时刻 ‖δX (τ0) ‖ 充分小, 于是1维的Lyapunov指数可定义为

在n维连续系统中, δX (τ) 在每个基底上有分量, 每一个分量均可按上式求出一个λ, 因此共存在n个Lyapunov指数λi, 称为Lyapunov指数谱。当任意选取矩阵δX (τ) 时, Lyapunov指数以概率1 可能取得最大值, 如果其中最大的Lyapunov指数λmax>0, 则该系统一定存在混沌运动。因此, 只要计算出系统的最大Lyapunov指数, 就可以判断系统是否处于混沌状态。

3 计算与分析

图4为柔性转子系统简图, 转轴由两个型号相同的滚子轴承支承, 转轴中有刚性圆盘, 该圆盘可施加不平衡力。

表1所示为滚子轴承的结构参数, 其中, 轴承的弹性模量为204GPa, 泊松比为0.3, 阻尼为200N·s/m, 轴承载荷为{0, 2000N, 2000N, 0, 0}。

柔性转子系统的结构参数如表2所示, 转轴被划分为4 个轴段, 转子系统划分为5 个节点。转轴的弹性模量为204GPa, 泊松比为0.3。

对图4所示的柔性转子系统进行计算, 判断圆盘节点处的混沌行为, 并分析轴承径向间隙、圆盘不平衡力、转轴刚度比等参数对系统混沌行为的影响规律。

3.1 径向间隙影响分析

假设系统不平衡力为0, 则轴承径向间隙对系统混沌特性的影响如图5所示。分别取40μm、60μm和80μm的径向间隙计算系统的最大Lya-punov指数。

系统最大Lyapunov指数小于0, 表明系统运动是稳定的;当最大Lyapunov指数等于0时, 系统为倍周期或拟周期分叉运动;当最大Lyapunov指数大于0时, 系统为混沌运动。

为了验证计算的准确性, 使用Runge-Kutta法进行非线性动力学方程组求解, 得到转子系统响应的频谱图和Poincaré截面。从图5可知, 当径向间隙为80μm、转速为4200r/min时, 系统的最大Lyapunov指数等于0 (图中fVC为滚子轴承的变刚度振动频率) 。图6所示为图5工况下圆盘节点处的频谱图和Poincaré截面, 图形表明系统为拟周期振动。

图5 表明, 径向间隙为80μm、转速为8000r/min时, 系统的最大Lyapunov指数大于0, 图7所示为对应该工况的圆盘节点处的频谱图和Poincaré截面, 其表明系统为混沌运动。

图5还表明, 当径向间隙为40μm时, 系统的最大Lyapunov指数小于0, 未出现混沌行为;当径向间隙为60μm时, 在0~2500r/min、6600~7200r/min转速下的最大Lyapunov指数大于0, 系统出现混沌行为;当径向间隙为80μm时, 在0~3100r/min、5000~5400r/min、7500~8500r/min转速下出现混沌行为。

3.2 不平衡力影响分析

图8所示为在其他条件不变前提下的不平衡力对系统混沌特性的影响。分别取圆盘的不平衡质径积为5×10-6kg·m、1×10-5kg·m、2×10-5kg·m, 计算系统的最大Lyapunov指数。

从图8 可知, 当不平衡质径积为5×10-6kg·m时, 在0~2000r/min, 系统出现混沌行为;当不平衡质径积为1×10-5kg·m时, 在0~3700r/min, 系统出现混沌行为;当不平衡质径积为2×10-5kg·m时, 在0~4400r/min, 系统出现混沌行为。可以看出, 不平衡力的存在会改变系统的混沌行为, 系统出现混沌的转速及范围随着不平衡力的增大逐渐增大。

3.3 刚度比影响分析

定义柔性转子系统的刚度比为

式中, Kshaft为柔性转轴的刚度;Kbearing为支承轴承刚度。

在其他参数不变的情况下, 分别取长度为500mm、400mm、300mm的转轴进行刚度比计算, 其对应所得的刚度比分别为7.38×10-3、1.43×10-2、3.31×10-2, 据此, 可得出刚度对于柔性转子系统动力学特性的影响。

图9所示为不同刚度比时圆盘节点处的振幅随转速变化曲线, 可以看出, 随着刚度比的增大, 振幅峰值对应的转速逐渐增大, 即临界转速增大;振幅峰值也随着刚度比增大而增大, 且临界转速附近的峰值逐渐变多, 可以看出随着转轴刚度增大, 轴承非线性振动对转子系统的影响增大, 导致系统响应复杂。

图10所示为不同刚度比时最大Lyapunov指数随转速的变化曲线。 当刚度比为7.38×10-3时, 在0~800r/min出现混沌行为。当刚度比为1.43×10-2时, 在0~1200r/min、5200~5400r/min出现混沌行为。当刚度比为3.31×10-2时, 在0~3200r/min、4800~5900r/min出现混沌行为。

可见, 随着转轴刚度增大即刚度比增大, 混沌运动的区间发生改变。随着刚度比的增大, 系统的混沌区间增大, 轴承引起的转子系统非线性行为明显。对比柔性转子和刚性转子的混沌特性, 发现轴承的非线性接触力对刚性转子系统的混沌特性影响大于柔性转子。

4 实验

为验证本文方法对柔性转子系统动力学行为预测的准确性, 使用图11所示的高速滚子轴承柔性转子实验器进行测试。动力装置的高速电主轴, 最高转速可达24 000r/min, 实验滚子轴承参数如表1所示, 实验过程中滚子轴承承受的径向载荷为2000N。采用非接触式的电涡流传感器测量滚子轴承及圆盘处的径向振动位移。

图12所示为转子圆盘处振动位移幅值计算结果与实验结果对比, 可以看出, 计算得到的一阶临界转速为4900r/min, 实验得到的一阶临界转速为4700r/min, 与计算结果较为接近;实验测试振幅与计算结果也较为接近, 证明了本文方法的准确性。

5 结论

(1) 轴承径向间隙是影响转子系统非线性振动特性的重要参数。随着轴承径向间隙的增大, 系统的混沌区间逐渐增大、变多。

(2) 不平衡力的存在对系统混沌行为也有较大的影响。系统出现混沌的转速及范围随着不平衡力的增大逐渐增大。

混沌行为的判据 篇4

电力系统负荷预测是电力系统调度、规划、供电等管理部门的基础工作,准确、有效的负荷预测不仅可以合理安排电网内部机组的启停、保持电网安全稳定地运行,还可以减少一些不必要的储备容量,合理安排检修计划,从而保证了正常的生产,有利于经济效益和社会效益的提高[1]。过去的几十年来,国内外学者将各种预测方法和模型运用到电力系统短期负荷预测中,使预测精度得到了很大的提高。文献[2]把粗糙集和神经网络结合建立短期负荷预测模型,采用粗糙集理论对各种影响负荷预测的因素变量进行识别,以此确定预测模型的输入变量;在此基础上通过属性约简和属性值约简获得推理规则集,再以这些推理规则构筑神经网络预测模型,并采用附加动量项的BP学习算法对网络进行优化,但是该方法没有对工作日和休息日的负荷预测加以区分,预测精度不够;文献[3]采用改进的粒子群算法和BP神经网络结合,提出了在算法迭代过程中,每个粒子会额外生成与迭代次数相同的粒子,并与当前粒子同方向不同速度飞行,利用适应度值保存粒子历史最优值。虽然也改善了粒子多样性,但这种方法是以显著增加计算量和牺牲系统内存为代价;文献[4]使用PSO算法优化基函数中心和宽度,再用最小二乘法确定隐含层与输出层间的权值,最后将改进算法应用于时间序列的预测中。但该方法初始粒子群随机产生,会导致算法收敛速度的不确定性,降低算法的平均收敛速度;文献[5]采用混沌神经网络对短期负荷进行预测,但是仅运用混沌时间序列分析作为神经网络日峰值预测模型选择最佳嵌入维数和延迟时间的必要理论依据,其不足之处是仅局限于日峰负荷预测,同时对于混沌网络的权值和阈值的确定较为困难且速度慢。

本文采用量子化粒子群算法不仅通过全同粒子系改善了初始种群的质量,而且通过对粒子的全局最优值与粒子的局部最优值的比较,限制粒子陷入局部最小搜索状态,提高粒子的局部搜索能力,节省了搜索时间,使粒子能够快速地搜索到最佳位置,从而增强了算法的局部寻优能力和收敛速度及计算精度;利用优化后的粒子群算法确定混沌神经网络的权值和阈值,克服混沌神经网络参数确定难度大、速度慢的缺点。本文在建立负荷预测模型的时候,考虑了休息日和工作日的日负荷不同的特点建立新的预测模型,提高了预测的精度。

2 基本粒子群算法及其改进

2.1 基本粒子群优化算法描述

在基本粒子群算法中,种群是由n个粒子组成的,粒子i的信息表示为d维向量[6]。位置用xi=(xi,1,xi,2,...,xi,d)(i=1,2,...,n)表示,速度为vi=(vi,1,vi,2,...,vi,d)(i=1,2,...,n),pi=(pi,1,pi,2,...,pi,d)表示第i个粒子的最优位置,其他向量类似。速度和位置更新公式为:

其中,vi,j(t)是粒子i在第t次迭代中第j维的速度;c1,c2是加速系数(或称学习因子),控制粒子群向全局最好粒子和个体最好粒子方向飞行的最大步长,适当的c1,c2取值能够加快粒子群的收敛速度并且使粒子群不易陷入局部最优;r1,r2是[0,1]之间的随机数;xi,j(t)是粒子i在第t次迭代中第j维的当前位置;pi,j是粒子i在第j维的个体极值点的位置;pg,j是种群在第j维的全局极点的位置。粒子的每一维速度v控制在(vmin,i,vmax,i)之间。vmax,i如果过大,粒子将会飞离最优解,太小将会陷入局部最优。假设将搜索空间的第d维定义为区间(xmin,i,xmax,i),每一维都用相同的设置方法。

2.2 量子行为粒子群优化算法

在基于量子行为的粒子群优化算法中,粒子的量子态通过波函数ψ(r,t)来表示。当ψ(r,t)确定后,粒子的所有力学分量和测值概率都可以确定,即xi,j(t)是由ψ(r,t)2决定的。在量子力学理论中,将属性相同的粒子称为全同粒子,由于全同粒子系具有交换对称性的特点,使得波函数具有很大的限制。一般来说,全同粒子系的波函数ψ(q1,q2,...,qn)不一定表示粒子pi,j的本征态,所有的pi,j处于完全平等的地位。然而,所有pi,j的共同本征态是存在的,即是完全对称波函数和完全反对称波函数。

全同粒子系的波函数约束条件为:

根据蒙特卡洛方法,粒子群的运动等式可以转化为式(4),另外引入式(5)和式(7)。

其中,C为常数因子;p,α, 分别根据式(5)、(6)和(8)求取;u∈[0,1],为随机数据;α为t变量收缩因子,随着时间的变化而变化。

其中,α1和α2为t变量收缩因子的初值和终值;Tmax表示最大迭代的次数。根据经验数据,通常取α1=2.5,α2=0.5,因此α∈(0.5,2.5)。pg,j表示每一个粒子全局搜索的最佳位置;pi,j表示每个粒子局部搜索到的最优位置;M表示种群的大小。

3 混沌神经网络模型

混沌神经元结构图如图1所示。

其中,vij、wij分别是第j个神经元的输入连接权值和反馈连接权值,xj(t+1)是神经元的输出[7]。

混沌神经元的输出函数为:

其中

式中,a表示神经元之间的联接强度,也称耦合因子。

在混沌神经网络结构中,包括输入层、隐层和输出层。其中,隐层的每一个神经元都会受到外部输入和内部反馈的影响,通过不停地调节神经元的权值和阈值,得到合适的混沌神经网络模型[8]。

混沌神经网络的输出函数为:

其中,xi是单个混沌神经元的输出值;wo是输出层神经元的阈值;wi是输出层神经元的权值;n1为隐层神经元的个数[9]。假设网络外部输入时间序列为u(t),隐层输出为o(t),网络输出为y(t),混沌网络表示为:

f1采用Sigmoid函数,即:y=1/[1+exp(-x)];f2采用线性函数:1W、2W和HW分别为输入层至隐层、隐层至输出层以及隐层节点之间的连接权矩阵。

4 基于量子行为粒子群优化算法-混沌神经网络负荷预测

4.1 基本原理

本文采用的混合算法中,将粒子群的位置向量x作为混沌神经网络的节点间连接权值和阈值,在每次的迭代过程中,利用优化后的粒子群算法求出权值和阈值,然后利用混沌网络,求出对应的权值和阈值的实际输出值fk(k=1,2,...,n)(n是神经网络输入输出的样本对数)。

粒子的适应度函数为:

式中,yk是混沌神经网络的目标输出;fk是混沌神经网络的实际输出。

4.2 混合算法模型

针对电力系统的负荷具有周期性的特征,同时工作日和休息日的日负荷不同的特点,本文采用的模型如图2所示,采用多输入、单输出[10]。

对于混沌网络的隐层节点数的确定采用经验公式:

其中,n1为输入层节点数;n2为输出层节点数;N为修正值。根据多次实验结果,同时保证运算的速度,当n=5时,运算速度和结果的误差能够满足需要。本文中,取n=5。

4.3 算法分析

在混合算法中,首先是将神经网络的权向量和阈值作为粒子群搜索空间中位置元素,然后应用粒子群优化算法计算出神经网络的权向量和阈值,即求出每一个粒子相应的实际输出值ok(k=1,2,...,n;n是神经网络输入输出的样本对数)。第i(i=1,2,...)个粒子的适应度函数为:

其中

其中,yk是神经网络的目标输出。

本文采用平均绝对百分误差EM和均方根误差ER作为评估指标[11]。

4.4 粒子群-混沌神经网络混合算法的流程

(1)根据网络的输入和输出关系,初始化混沌网络的拓扑结构[12]。确定粒子的初始位置xi,j(0)和速度,确定粒子数M、最大允许迭代次数Tmax、加速系数c1和c2;

(2)如果是基本粒子群优化算法则用式(1)和式(2)对每一个粒子的速度和位置进行更新;如果是具备量子行为的粒子,采用改进粒子群优化算法,用式(7)和式(8)分别确定每个粒子的全局最优位置、局部搜索位置;

(3)根据优化后的粒子群算法,求出混沌神经网络的权值和阈值;

(4)利用混沌神经网络计算出每个粒子对应的个体极值,将粒子群中个体极值最好的作为全局极值。记录该粒子的序号,用gbest(全局极值点)表示最好粒子的当前位置;

(5)根据粒子的适应度函数,计算每一个粒子的适应度值。如果粒子的适应度值优于该粒子的个体极值,则将pbest(个体极值点)设置为该粒子的位置,同时对粒子的个体极值进行更新。当全部粒子的个体极值优于此时的全局极值时,将gbest设置为该粒子的位置,记录该粒子的序号,同时对全局极值进行更新;

(6)判断是否满足流程结束条件。如果当前位置满足预定要求(迭代次数达到了给定的最大次数或达到最小误差要求)时,则停止迭代,输出最优解;如果不能满足结束条件,转到步骤(2)。

4.5 数据的归一化处理

为了确保输入量具有较好的作用,选用Sigmoid函数中间段的函数关系,从而避开其两端的饱和区域,必须对神经网络的输入量进行归一化处理[13]。

t时刻负荷数据采用如下归一化公式:

在输出层则用式(22)重新换算回负荷值:

式中,Lmax和Lmin分别为训练样本集中负荷的最大值和最小值。

5 应用实例及结果

本文预测模型中混沌神经网络的反馈过程是通过循环实现的,其停止的条件用精度来判断,即如果A(t)-A(t-1)

本文结合某地的实际情况,对其某日24h整点的电力负荷分别采用量子粒子群算法、混沌学习算法和本文提出的量子粒子群优化-混沌神经网络算法进行预测,评估指标对比情况如表1所示,负荷预测结果如表2所示。

由表1可知,量子粒子群-混沌学习算法在训练550次左右的EM值已经小于量子粒子群算法和混沌学习算法训练2500次的EM值。量子粒子群-混沌神经网络算法在训练1500次的EM值也小于量子粒子群算法和混沌学习算法训练5000次时的EM值,所用时间前者80s,后两者的时间分别为213s和256s。可见在收敛性和训练速度上,本文采用的混合算法优势明显。三种算法预测结果与实际值的平均百分绝对误差对比图如图4所示。

对比表2中预测结果和相对误差可知,采用量子粒子群-混沌神经网络算法训练400次时的负荷预测结果精度已经好于采用量子粒子群算法(3000次)和混沌学习算法(2000次)时的预测精度,表明本文采用的预测方法和模型在预测精度和速度方面,明显好于以上两种算法。

从图5可以看出,本文采用的量子粒子群-混沌神经网络的混合算法预测结果相对误差控制在4%以内,且误差波动较小。预测精度比量子粒子群算法和混沌学习算法要好很多。

从图6中可以看出,当迭代次数达到900时,量子粒子群-混沌神经网络算法的适应度函数就基本达到稳定。而量子粒子群算法和混沌学习算法迭代次数达到1700和2000次左右时候才达到稳定。将粒子群的适应度函数设定为训练误差,适应度函数越大,输出误差越大。由此可见,本文采用的混合算法模型的辨识精度远高于其他两种算法模型的辨识精度,表明本文采用的模型更加实用。

6 结论