System size resonance in coupled noisy systems and in the Ising model
- 格式:pdf
- 大小:417.63 KB
- 文档页数:4
第27卷㊀第12期2023年12月㊀电㊀机㊀与㊀控㊀制㊀学㊀报Electri c ㊀Machines ㊀and ㊀Control㊀Vol.27No.12Dec.2023㊀㊀㊀㊀㊀㊀DAB 级联单相逆变器系统的阻抗特性及稳定性分析刘欣,㊀袁静,㊀高鑫波(华北电力大学电气与电子工程学院,河北保定071003)摘㊀要:针对双有源桥(DAB )直流变换器级联单相并网逆变器系统因阻抗失配而造成系统发生振荡失稳的问题,通过建立DAB 和单相逆变器的输出和输入阻抗模型,基于阻抗分析法对级联系统的交互稳定性进行了分析㊂首先,推导采用双环控制策略的前级DAB 输出阻抗模型和考虑锁相环影响的后级逆变器直流侧输入阻抗模型,并通过扫频法验证其准确性㊂在此基础上,建立二者阻抗交互模型,详细分析了DAB 反馈控制器的PI 参数对其输出阻抗频率特性和级联系统稳定性的影响,并据此提出一种DAB 控制参数的优化设计方法,在兼顾动态性能的同时提升了级联系统的稳定性㊂最后,通过仿真算例验证了阻抗模型的准确性,分析了结论的正确性以及稳定性改善方法的有效性㊂关键词:级联系统;稳定性;阻抗重塑;双有源桥;单相并网逆变器;阻抗模型DOI :10.15938/j.emc.2023.12.001中图分类号:TM46文献标志码:A文章编号:1007-449X(2023)12-0001-11㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀收稿日期:2023-03-17作者简介:刘㊀欣(1980 ),男,博士,副教授,研究方向为新能源发电系统建模与控制㊁电力电子系统电磁兼容和瞬态特性;袁㊀静(1997 ),女,硕士研究生,研究方向为电力电子变流器建模与控制;高鑫波(1999 ),男,硕士研究生,研究方向为电力电子变流器建模与控制㊂通信作者:袁㊀静Impedance characteristics and stability analysis of DAB cascadesingle-phase inverter systemLIU Xin,㊀YUAN Jing,㊀GAO Xinbo(School of Electrical and Electronic Engineering,North China Electric Power University,Baoding 071003,China)Abstract :Aiming at the problem of oscillation instability of the dual active bridge (DAB)DC-DC con-verter cascaded single-phase grid-connected inverter systems due to its impedance mismatching,the out-put and input impedance models of DAB and single-phase inverter were established,and the interaction stability of the cascade system was analyzed based on impedance analysis method.Firstly,the output im-pedance model of the front-stage DAB using the double-loop control strategy and the DC-side input im-pedance model of the back-stage inverter considering the influence of the phase-locked loop were derived,and the accuracy of the models were verified by frequency sweep method.Based on this,the impedance interaction model between the two was established.Additionally,the effects of PI parameters of DAB feedback controller on its output impedance frequency characteristics and cascade system stability wereanalyzed in detail,and the optimal design method of DAB control parameters was proposed accordingly,which improves the stability of the cascade system while taking into account the dynamic performance.Fi-nally,the simulation examples verify accuracy of the impedance model,correctness of the analytical con-clusions and effectiveness of the stability improvement method.Keywords :cascaded system;stability;impedance reshaping;dual active bridges;single-phase grid-con-nected inverters;impedance model0㊀引㊀言在光伏系统㊁蓄电池㊁超级电容,车网互联(ve-hicle to grid,V2G)等交流并网型储能系统中,通常需要使用两级式DC /AC 变换器实现并入交流电网和双向功率控制的功能[1]㊂其中,双有源桥变换器由于具有高功率密度㊁电流隔离㊁能量双向传输和易实现零电压开关等优点[2-4],很好地适应了交流并网型储能系统的需求,是第一级DC /DC 变换器的理想选择,而单相逆变器用于与电网连接㊂基于双有源桥(dual active bridge,DAB)变换器的两级式DC /AC 变换器的典型电路拓扑如图1所示㊂该拓扑整体结构简单,易于实现,控制方法较为成熟,被大量应用于电动汽车充电桩领域[5-8]㊂然而,由于变换器复杂的输入输出特性以及级联结构的存在,尽管两级变换器在单独运行时能保持稳定,但子系统之间的相互作用可能会使系统性能下降,导致直流母线产生电压振荡,以至于系统崩溃[9]㊂因此,通过稳定性分析㊁合理参数调整㊁控制优化等方法改善级联系统的稳定性和可靠性是当今研究的一个热点与难点问题[10-12]㊂图1㊀两级式DC /AC 变换器主电路拓扑及控制框图Fig.1㊀Main circuit topology and control block diagram of two-stage DC /AC converter㊀㊀基于阻抗的Nyquist 阻抗匹配原则[13]已经被广泛应用于各类级联系统的交互稳定性的研究中㊂准确的阻抗模型对于级联系统稳定性分析是必要的㊂目前,常用的逆变器阻抗建模方法包括谐波线性化法[14-16]和dq 坐标系下的阻抗建模法[17-18]㊂谐波线性化将系统视为2个单输入单输出系统,主要用于分析三相系统的谐波稳定性;而dq 阻抗建模法通常将电气量转变为d 轴和q 轴分量,以便单独控制有功和无功功率,有利于在稳态工作点处进行小信号分析㊂文献[19]在dq 坐标系下推导了使用不同控制策略的三相并网逆变器的直流侧输入阻抗模型,此方法适用性较强,但并未应用到单相逆变器系统中㊂文献[20]提出一种基于二阶广义积分器(second order generalized integrator,SOGI)的dq 坐标系下单相整流器的阻抗建模方法,但此方法并未推广到单相并网逆变器的阻抗建模中㊂由于阻抗相互作用是造成两级式DC /AC 级联系统失去稳定的根本原因,可以通过重塑源变换器或者负载变换器的总线端口阻抗来提高系统的稳定性㊂为了达到这一目的,学者们提出了多种方法,包括无源阻尼法[21-23]和有源阻尼法[24-26]㊂其中,无源阻尼法需要引入附加无源元件,以改变变换器的阻抗特性,但附加阻尼电路会增加硬件成本,降低变换器效率;有源补偿法具有成本低㊁不增加损耗的优点,因而被广泛用于基于DAB 变换器的级联系统阻抗匹配优化设计中㊂文献[27]采用有源阻尼的优化思路对LC -DAB 级联系统进行阻抗重塑,提出基于一次侧电容电压的并联虚拟阻抗和一次电流串联虚拟阻抗控制策略,从而使得级联系统在全功率范围内均能稳定运行;文献[28]研究基于DAB 的储能系统稳定性,提出在窄带范围内对负载变换器DAB 的输入阻抗进行重塑,在提高稳定性的同时保证系统动态性能良好;文献[7]研究了用于电动汽车双向充放电的DAB 级联单相并网电压源变换器(voltage source converter,VSC)系统的阻抗稳定性,提出一种基于虚拟电阻的有源阻尼方法,以改变2电㊀机㊀与㊀控㊀制㊀学㊀报㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀第27卷㊀VSC的输出阻抗,提高级联系统在各种工作模式下的稳定性;文献[29]面向DAB级联三相VSG系统,通过构建与DAB转换器的输入阻抗并联或串联的虚拟阻抗以增加DAB输入阻抗幅值,从而满足稳定性准则㊂文献[30]针对具有电压调整单元的DAB 变换器提出一种基于超前-滞后的阻抗优化调节器用以抑制输出阻抗谐振尖峰,提升了系统运行可靠性,并优化了电流应力㊂总的来说,上述级联系统的稳定性增强方法都需要增加附加的控制过程,从而不可避免地增加了模型的复杂度,其设计方法仍存在进一步简化的空间㊂而DAB变换器的输入阻抗会受到其反馈控制器的影响,揭示二者之间的关联有助于简化阻抗匹配优化设计,但此方面的相关研究较少,并且缺乏深入的理论分析㊂针对上述问题,本文对双有源桥DC/DC变换器与单相并网逆变器组成的级联系统进行阻抗建模并进行稳定性分析㊂首先,建立采用双环控制策略的DAB输出阻抗模型和采用解耦电流控制策略的单相并网逆变器直流端输入阻抗模型,并将锁相环的相位波动考虑在内,通过扫频法验证阻抗模型的正确性㊂随后,建立阻抗交互模型,从理论上分析DAB变换器的PI参数对其输出阻抗波形的影响,结合Nyquist图和闭环根轨迹进一步讨论关键参数与系统稳定性之间的关联㊂分析结论表明,调节DAB电压外环比例系数可直接调节级联系统稳定性,基于此,提出通过优化DAB变换器的电压外环比例系数提高级联系统稳定性的方法,该方法无需任何额外的补偿器或控制回路,在兼顾系统动态性能的同时,有效实现了基于DAB的交直流级联系统的稳定性增强㊂MATLAB/Simulink仿真算例验证了稳定性改善方法的有效性㊂1㊀级联系统阻抗建模变换器阻抗的精确建模是稳定性分析的基础㊂图1所示的控制框图为级联系统的常规控制方案,其中,DAB变换器负责控制直流母线电压的稳定,单相并网逆变器负责控制功率输出[31-32]㊂本节将分别给出DAD输出阻抗和单相并网逆变器的直流侧输入阻抗的建模过程㊂1.1㊀DAB变换器输出阻抗建模DAB变换器的拓扑及控制方案如图1中左面虚线框所示㊂其输出功率[33-34]可表示为P=nV in v busL o f s dϕ(1-2dϕ)=v bus i2⓪㊂(1)式中:n为变压器变比;V in为DAB输入电压;v bus为输出电压;L o为变压器等效电感;f s为开关频率;dϕ为变压器两侧H桥输出电压之间的相移量(dϕ=ϕ/2π);i2为副边H桥输出电流, i2⓪表示其平均值㊂经小信号分析可得i2与占空比dϕ的关系为G i2d=i^2d^ϕ=nV in Lo f s(1-4Dϕ)㊂(2)式中符号^表示变量的小信号形式㊂采用内环电流加外环电压的双环控制模式㊂将控制器的内环传函记作G c1(s),外环传递函数记作G c2(s),其中:G c1(s)=k pi+k ii s;G c2(s)=k pv+k iv s㊂将负载变换器阻抗等效为R,则DAB控制回路小信号模型如图2所示,图中LPF为一阶低通滤波器,用于实现20dB/dec的环路增益[35](H LPF(s)= 1/(s/ωLPF+1),其中ωLPF为低通滤波器的截止频率)㊂图2㊀DAB控制回路小信号模型Fig.2㊀Small signal model of DAB control loop根据上述控制框图,得到DAB的输出阻抗为Z out_DAB=v^bus-i^bus=1C bus s+G c1G x㊂(3)式中G x=G c2G i2d1+G c2G i2d H LPF㊂1.2㊀单相并网逆变器直流侧输入阻抗建模基于SOGI的锁相环(PLL)模型如图3(a)所示㊂图中,v为自公共耦合点(PCC)电压(将其本身视为静止坐标系下的α轴分量,β轴虚拟分量与之垂直)㊂SOGI的传递函数为H e(s)=K SOGIω1ss2+K SOGIω1s+ω21㊂(4)式中:ω1为电网工频;K SOGI为闭环系数㊂在小扰动下,PLL输出与PCC实际相位存在相位差Δθ,其将导致控制系统中的各变量与功率系统中的相应变量存在差异㊂为以示区分,文中带有上标s的变量表示 电气量 ,带有上标c的变量表示 控制量 ㊂为了简化表达式,将成对变量以矢量形式编写,例如v s dq表示[v s d v s q]T,另外,变量的大写符3第12期刘㊀欣等:DAB级联单相逆变器系统的阻抗特性及稳定性分析号表示其自身静态工作点㊂图3㊀基于SOGI 的PLL 模型Fig.3㊀SOGI-based PLL model根据图1可得系统功率方程为:(Z L +Z g )i ^s dq =D dq v ^bus +d ^sdq V bus ;i ^bus=12(D T dq i ^s dq +I T dq d ^sdq )㊂}(5)式中:i ^s dq =[i ^s d i ^s q ]T 和d ^s dq =[d ^s d d ^s q ]T分别为交流侧电流与占空比的dq 轴分量构成的列向量;Z L =sL f +R f -ωL f ωL f sL f +R f éëêêùûúú;Z g =sL g +R g -ωL g ωL g sL g +R g éëêêùûúú;L f 和R f 为滤波电感及其等效内阻;L g 和R g 为电网内阻抗;i ^bus 为逆变器直流侧输入电流㊂将图3中Park 变换框图T θ1前移,得到其等效控制框图如图3(b)所示,图中:H edq (s )=A B-B A[];A =[H e (s +j ω1)+H e (s -j ω1)]/2;B =[j H e (s +j ω1)-j H e (s -j ω1)]/2㊂根据图3(b)可推导PCC 电压的 控制量v ^c dq与 电气量v ^sdq之间的关系为v ^c dq =G v PLL v ^sdq ㊂(6)式中:Gv PLL=v ^c dq v^s dq=A -V sqG s B B +V sqG s A -B +V sd G s B A -V sd G s Aéëêêùûúú;G s 为PLL 输出角度与PCC 电压q 轴分量的关系式;G s =sk p_PLL +k i_PLLs +V s d (sk p_PLL +k i_PLL ),k p_PLL 和k i_PLL 为锁相环PLL 的PI 参数㊂同理可得输出电流与占空比的 控制量 与 电气量 的小信号关系为:d ^sdq=d^cdq+G dPLL v ^s dq;i^c dq=Hedq i^s dq+GiPLL v ^s dq㊂}(7)式中:G d PLL =D s qG s B -D s qG s A -D sd G s B D sd G s Aéëêêùûúú;Gi PLL=-I sq G s B I s q G s A I sd G s B-I sd G s Aéëêêùûúú㊂令:H i =k p_INV +k i_INV /sk p_INV+k i_INV /s éëêêùûúú,其中:k p_INV 和k i_INV 为逆变器电流控制器的PI 参数;G ci=k p_INV +k i_INV /s ωL f-ωL fk p_INV+k i_INV /s éëêêùûúú㊂将解耦电流控制策略与PCC 电压前馈结合,得到考虑锁相环影响的逆变器控制回路的小信号模型如图4所示㊂图4㊀PLL 影响下电流控制回路小信号模型Fig.4㊀Small-signal model for current control loopwith PLL根据图4,得到逆变器控制部分的方程为d ^s dq =[(G v PLL -G ci G i PLL +V bus G dPLL )Z g -G ci H edq ]i ^s dq /V bus ㊂(8)联立式(5)㊁式(8)可得单相并网逆变器直流侧输入导纳为Y in_INV =i ^busv ^bus=12V busI T dq (Z L +Z g )+12D T dq[]㊃([Z L +G ci H edq -G PLL_V Z g ]-1D dq )-12V bus I Tdq D dq㊂(9)式中G PLL_V =G v PLL -G ci G i PLL +V bus G dPLL -E ,其中E为单位矩阵㊂相应的单相并网逆变器直流侧输入阻抗为Z in_INV =1/Y in_INV ㊂(10)1.3㊀阻抗模型的仿真验证基于MATLAB /Simulink 平台搭建了DAB 级联单相并网逆变器的仿真模型,采用扫频法对2个级联子系统的输出和输入阻抗模型分别进行验证,仿真参数如表1所示㊂图5给出了仿真扫频与理论模型的对比结果㊂可以看出,在1~10000Hz 频段,所得阻抗模型与扫频结果吻合较好,验证了所推得的DAB 输出阻抗和单相逆变器输入阻抗模型的正确性㊂4电㊀机㊀与㊀控㊀制㊀学㊀报㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀第27卷㊀图5㊀级联系统阻抗模型Fig.5㊀Impedance model of cascade system表1㊀级联系统电路参数Table1㊀Parameters of cascade system㊀㊀㊀参数数值DAB直流侧输入电压稳态值V in/V400直流母线电压稳态值V bus/V400直流母线电容C bus/μF1500 DAB变压器等效电感L o/μH30变压器变比n1ʒ1 DAB开关频率f s/kHz20低通滤波器截至频率ωLPF/(rad/s)4000π逆变器并网电压有效值V g/V220逆变器输出功率稳态值P/kW10逆变器滤波电感L f/mH及等效内阻R f/mΩ10,50电网内电感L g/mH及内电阻R g/mΩ1,152㊀级联系统稳定性分析2.1㊀级联系统阻抗交互模型级联系统的稳定性不仅取决于变换器各自的稳定性,还决定于源变换器(本文为DAB)输出阻抗与负载变换器(本文为单相逆变器)输入阻抗二者交互作用的影响㊂将DAB视为电压源,逆变器视为电流源,二者构成的级联系统阻抗相互作用示意图如图6所示㊂图6㊀级联系统等效阻抗示意图Fig.6㊀Equivalent impedance diagram of cascade system根据图6,可知级联系统开环传递函数为T m=Z out_DABZ in_INV㊂(11)式中T m也称为系统小环路增益㊂根据Middlebrook 判据[13],当源变换器和负载变换器各自稳定,并且系统的小环路增益T m满足Nyquist稳定判据时,该级联系统方是稳定的㊂图7为DAB输出阻抗和逆变器输入阻抗伯德图㊂由于逆变器采用恒功率控制,因此,除50Hz频点外,在f<f c3(f c3为逆变器电流控制器的截止频率)频率范围内逆变器直流端输入阻抗呈现阻值为-V2bus/P的负电阻特性;在f>f c3频率范围内呈现电感性质㊂而50Hz频点是一个特殊点,其阻抗幅值几乎为0,相位跃变到0㊂虽然逆变器与DAB的阻抗在50Hz频点处容易产生交叉,但二者相位之差小于180ʎ,因此不影响系统稳定性㊂此外,DAB输出阻抗在f<f r频段(f r为DAB输出阻抗谐振频率)呈现电感特性,在f>f r频段呈现电容特性㊂这使得DAB输出阻抗具有类似LC滤波器的阻抗特性㊂图7㊀级联系统阻抗伯德图Fig.7㊀Impedance Bode diagram of cascade system综合以上阻抗特性可知,DAB输出阻抗的谐振峰以及逆变器在低频段的负阻抗特性是导致交直流级联系统稳定性降低的主要原因㊂一旦DAB输出阻抗的谐振峰与逆变器输入阻抗发生交叉,就会因5第12期刘㊀欣等:DAB级联单相逆变器系统的阻抗特性及稳定性分析相位裕度无法满足稳定条件而造成系统振荡失稳㊂2.2㊀DAB 变换器的控制参数分析由图7可知,平抑DAB 输出阻抗的谐振峰将有效提高级联系统稳定性㊂为了达成这一目的,本节将详细分析DAB 反馈控制器的PI 参数与谐振峰之间的关联,为级联系统的稳定性分析及控制器参数优化设计奠定基础㊂当DAB 电流内环截止频率与一阶低通滤波器LPF 带宽相等时,经控制器定量设计可得电流控制器比例系数k pi 为0㊂将k pi =0代入式(3),并且忽略含有C bus 和T LPF 的高阶项,整理得到DAB 输出阻抗的简化表达式为Zᶄout_DABʈ1G i2d k ii s (s +G i2d k ii )C bus s 2+k pv s +k iv㊂(12)图8给出了DAB 输出阻抗的理论模型和简化模型的对比图㊂可以看出,在1~200Hz 频率范围内,二者阻抗模型基本吻合,结合图7可知,影响系统稳定性的频段为几十赫兹,因此说明上述简化模型可胜任稳定性分析需求㊂图8㊀DAB 理论模型和简化模型对比Fig.8㊀Comparison of theoretical and simplified Bodediagrams of DAB设定DAB 电流内环截至频率f c1为2000Hz,相位裕度P m1为45ʎ,同时电压外环截至频率f c2为20Hz,相位裕度P m2为45ʎ时,经设计所得DAB 的控制参数如表2所示㊂表2㊀DAB 控制器PI 参数Table 2㊀PI parameters of DAB controller㊀㊀㊀㊀参数数值电流控制器比例系数k pi 0电流控制器积分系数k ii 30.443电压控制器比例系数k pv 0.102电压控制器积分系数k iv24.35㊀㊀将s =j ω代入式(12),得到DAB 阻抗的模值为|Z ᶄout_DAB (j ω)|=ωaω2(C bus ω2-k iv +ak pv )2+(ω2k pv -C bus ω2a +ak iv )2(-C bus ω2+k iv )2+ω2k 2pv㊂(13)式中a =G i2d k ii ㊂令Z ᶄout_DAB (j ω)虚部为0,得到谐振点频率为ω=G i2d k ii k ivC bus G i2d k ii -k pv㊂(14)根据式(13)和式(14)可得DAB 输出阻抗的谐振频率及谐振峰值分别与控制参数的关系曲线如图9所示㊂结合式(27)㊁式(28)和图9,可得如下结论:当电压外环比例系数k pv 增大时,谐振频率几乎不变,谐振峰值陡然降低;当电压外环积分系数k iv 增大时,谐振频率增大,谐振峰值维持不变;当电流环积分系数k ii 改变时,二者均基本不发生改变㊂上述分析表明,参数k pv 是平抑DAB 输出阻抗谐振峰的关键参数,而参数k iv 是改变谐振频点的关键参数㊂图9㊀谐振频率及谐振峰值与DAB 控制参数的关系曲线Fig.9㊀Relationship curves of resonant frequency andresonant peak with DAB control parameters为了佐证此结论,图10给出了不同控制参数下的DAB 输出阻抗伯德图㊂可以看出,当比例系数k pv 从0.02逐渐增大到0.4,且其余参数与表1和表2保持一致时,DAB 输出阻抗谐振峰值急剧减小,但谐振频点基本保持不变;当积分系数k iv 从10增大到120,且其余参数与表1和表2保持一致时,DAB 输出阻抗谐振频率逐渐增大,而谐振峰值几乎不变㊂6电㊀机㊀与㊀控㊀制㊀学㊀报㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀第27卷㊀图10㊀不同控制参数作用下DAB输出阻抗伯德图Fig.10㊀DAB output impedance Bode diagram of differ-ent control parameters综上所述,DAB电压外环控制器参数直接决定了其输出阻抗谐振峰值的大小及位置㊂其中,参数k pv与谐振峰幅值大小具有强相关性,适度增大参数k pv将显著降低DAB输出阻抗谐振峰,从而避免与逆变器输入阻抗发生交叉㊂据此可推断,参数k pv是作为影响级联系统稳定性的关键参数,对其进行优化设计可实现系统稳定控制,且设计过程也最为简单,相关分析及验证将在2.3节给出㊂2.3㊀DAB电压外环比例系数对系统稳定性的影响本节进一步讨论k pv对级联系统稳定性的影响㊂根据式(11)可知系统的特征方程为1+T m=0㊂(15)将式(3)和式(10)代入式(15),可得sC bus s2+(k pv s+k iv)G x+Z in_INV=0㊂(16)由于参数k pv直接体现在系统特征方程中,因此可结合基于闭环传递函数的根轨迹和开环传递函数的Nyquist图进行分析㊂对式(16)进行等效变换,保证特征方程不变,得到系统等效的开环传递函数为D(s)=k pv sG x Z in_INVs+Z in_INV(Cs2+k iv G x)㊂(17)根据式(17),得到当参数k pv从0逐渐变化至+ɕ时系统闭环传递函数的特征根在复平面的变化轨迹如图11所示㊂此时DAB电流内环控制参数与表2中相同,电压外环积分系数为98.3㊂可以看出,当k pv<0.0634时,级联系统存在右半平面极点,系统处于不稳定状态;当k pv>0.0634时,系统方可稳定;当k pv=0.0634时,复平面上出现位于虚轴上的闭环极点(0,ʃj251),说明系统处于临界稳定状态,这意味着系统中将会出现251rad/s(约40Hz)的振荡频率㊂图11㊀系统关于参数的k pv的根轨迹图Fig.11㊀Root trajectory diagram of the system with re-spect to the parameter k pv图12给出了此临界稳定状态下系统开环传递函数T m的Nyquist图,在此参数状态下,Nyquist曲线恰好穿越(-1,j0)点㊂分析结果说明,增大DAB电压外环比例系数k pv有助于增强级联系统稳定性,反之,将使级联系统稳定性变差㊂图12㊀系统开环传递函数的Nyquist图Fig.12㊀Nyquist diagram of the open-loop transferfunction为了验证上述分析结论,在MATLAB/Simulink 中搭建DAB与单相并网逆变器级联系统的仿真模7第12期刘㊀欣等:DAB级联单相逆变器系统的阻抗特性及稳定性分析型㊂电路参数如表1所示㊂图13给出了当其余参数保持不变,DAB 电压外环比例系数k pv 分别为2㊁0.258㊁0.0634和0.03时直流母线电压和交流侧输出电流的时域仿真波形㊂可以看出,当k pv 为2和0.258时,系统运行在稳定状态;当k pv 为0.0634时,系统处于临界稳定状态;当k pv 减小到0.03时,系统振荡失稳㊂这与图11中的参数根轨迹分析结果相符㊂图13㊀k pv 减小时直流母线电压和交流电流时域波形Fig.13㊀Waveforms of DC bus voltage and AC currentwhen k pv decreases取时间窗为0.2s,对图13中各个时间段的直流母线电压的时域波形进行频谱分析,所得结果如图14所示㊂可以看出,当k pv >0.0634时,直流母线电压主要含有直流分量和单相交直流系统中固有的二倍频分量;当k pv =0.0634时,在直流母线电压中出现可观的40Hz 频率分量,与图11中临界稳定状态下的系统振荡频率基本吻合;当k pv <0.0634时,直流母线电压中谐波分量杂乱繁多,系统失去稳定性㊂此外,还需特别说明的是,系数k pv 在影响系统稳定性的同时,也会影响系统动态响应速度㊂观察图11中根轨迹局部放大图可知,当k pv 小于0.258时,随着k pv 增大,主导极点的根轨迹从右半平面逐渐变化到左半平面并远离虚轴;当k pv 大于0.258时,根轨迹方向转变并逐渐靠近虚轴㊂因此,当k pv =0.258时,系统具有最佳的动态性能㊂若k pv 持续增大,越过最佳运行点,虽仍可保证稳定,但系统稳定速度将滞缓,这说明需兼顾稳定性和动态性能进行k pv 的参数设计㊂为了验证这一结论,图15给出了k pv 取值变化时系统有功功率波形的变化,从有功功率角度说明系数k pv 对系统稳定性及动态响应速度的影响㊂比较图15(a)㊁(b)和(c)可知,当0.0634<k pv <0.258时,系统动态响应速度随着k pv 的增大而加快,并且k pv 越大,系统稳定速度越快㊂比较图15(c)和图15(d)可知,当k pv 取2时,系统的稳定速度相较于图15(c)变慢,说明此时k pv 取值已越过了最佳运行点,进而验证了前述理论分析的正确性㊂图14㊀直流母线电压FFT 分析Fig14㊀FFT analysis of DC busvoltage图15㊀不同k pv 作用下的有功功率曲线Fig.15㊀Active power waveforms with different k pv8电㊀机㊀与㊀控㊀制㊀学㊀报㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀㊀第27卷㊀3㊀级联系统稳定性改善方法第2.3节中的分析结论表明,增大DAB 电压外环比例系数k pv 将显著提高级联系统稳定性㊂因此,当级联系统面临振荡失稳问题时,一种简单而可靠并且无需任何额外的补偿器或控制回路的稳定性改善方法为:增大DAB 电压外环比例系数k pv ㊂根据2.1节的分析,若要使系统满足稳定性要求,应保证增大k pv 后,DAB 输出阻抗的峰值小于V 2bus /P ,从而避免与逆变器输入阻抗发生交叉,并保证系统具有足够的相位裕度㊂本节将结合具体的仿真算例验证此稳定性改善方法的有效性㊂仿真算例中基本电路参数如表1所示,DAB 控制器的相关参数列于表2之中㊂图16给出了算例中直流母线电压和交流电流时域波形,图17则给出了与之对应的级联系统阻抗伯德图㊂如图16中0.2~0.5s 时间窗内波形所示,当系统传输功率为5kW 时,系统稳定运行,直流母线电压包含400V 的直流分量和二倍频分量㊂若传输功率增加为10kW,系统将发生振荡失稳,如图16中0.5~0.7s 的波形所示㊂图17中曲线Z in_INV1与Z in_INV2分别为功率改变前后逆变器输入阻抗伯德图㊂可以看出,负载的加重造成逆变器在低频段的阻抗幅值减小,因此与DAB 输出阻抗发生交叉㊂图16㊀稳定性改善前后的时域仿真波形Fig.16㊀Time domain simulation waveforms before andafter stability improvement为了改善系统稳定性,应当增大DAB 电压外环比例系数㊂图18为传输功率为10kW 时,级联系统关于参数k pv 的根轨迹曲线㊂可以看出,要想保证系统稳定运行,k pv 的取值必须大于0.0645,并且当k pv 取0.452时,系统具有最佳动态性能㊂观察图16中0.7~1.1s 时域波形可知,在0.8s 时,将DAB 电压外环比例系数调整为最佳参数0.452,其余参数保持不变,由于DAB 输出阻抗的谐振峰值降低,系统又重新恢复至稳定运行状态㊂图17㊀稳定性改善前后的级联子系统阻抗伯德图Fig.17㊀Impedance Bode diagram before and after sta-bility improvement of the cascadesubsystem图18㊀传输功率为10kW 时系统闭环根轨迹Fig.18㊀Closed-loop root trajectory of the system at10kW transmission power上述仿真算例进一步验证了稳定性改善方法的可行性㊂在系统控制器设计中,应当根据DAB 和单相并网逆变器的阻抗特性,利用阻抗伯德图和系统关于参数k pv 的闭环根轨迹进行直观判断,综合考虑系统的稳定性和动态响应速度,以确定适合的控制参数㊂4㊀结㊀论本文分别建立了DAB 输出阻抗模型和考虑锁相环相位波动影响的单相并网逆变器的直流端输入阻抗模型,提高了模型的准确度,并通过扫频法对阻抗模型进行验证;此外,通过理论分析获得了DAB 输出阻抗谐振频率及谐振峰值的计算公式,从原理9第12期刘㊀欣等:DAB 级联单相逆变器系统的阻抗特性及稳定性分析。
第51卷第10期电力系统保护与控制Vol.51 No.10 2023年5月16日Power System Protection and Control May 16, 2023 DOI: 10.19783/ki.pspc.221314考虑谐波耦合的多变流器并网系统建模及交直流谐波交互特性分析林顺富1,李 寅1,戴烨敏2,杨 帆1,边晓燕1,李东东1(1.上海电力大学电气工程学院,上海 200090;2.国网上海市电力公司青浦供电公司,上海 201799)摘要:新型电力系统中高比例可再生能源和高比例电力电子设备接入的特征带来异于传统电网下的谐波交互问题,多变流器拓扑的谐波耦合交互特性亟待研究。
首先,基于谐波状态空间(harmonic state space, HSS)对多变流器并网系统(multiple grid-connected-converter system, MGCCS)建立考虑谐波耦合的谐波传递函数矩阵模型,综合考虑了系统各控制环节对状态变量的影响以及变流器的级联、并联。
其次,基于所建HSS模型明确定义谐波耦合系数,并用于揭示多变流器拓扑的谐波耦合机理,分析级联、并联变流器谐波交互特性。
然后,应用谐波耦合系数量化分析MGCCS中滤波电感、电流环、锁相环等关键参数对系统谐波交互的影响。
最后,将HSS模型和Matlab/Simulink 模型、RT-LAB模型的结果进行对比,验证了所建HSS模型的精确性,以及谐波耦合系数理论应用于系统交直流谐波交互分析的有效性。
关键词:多变流器并网系统;交直流谐波交互;谐波状态空间;谐波耦合系数;并网电流Modeling of a multiple grid-connected-converter system considering harmonic coupling andanalysis of AC/DC harmonic interaction characteristicsLIN Shunfu1, LI Yin1, DAI Yemin2, YANG Fan1, BIAN Xiaoyan1, LI Dongdong1(1. College of Electrical Engineering, Shanghai University of Electric Power, Shanghai 200090, China; 2. Qingpu PowerSupply Company of State Grid Shanghai Municipal Electric Power Company, Shanghai 201799, China) Abstract: The characteristics of high penetration of renewable energy and high proportion of power electronic equipment in new power systems bring a harmonic interaction problem. This is different from traditional power systems. The harmonic coupling interaction characteristics of the multi-converter topology need to be studied. First, based on the harmonic state space (HSS), the harmonic transfer function matrix model considering harmonic coupling is established for the multiple grid-connected-converter system (MGCCS). The model comprehensively considers the influence of each control link of the system on the state variables as well as the factors of cascade and parallel converters. Second, the harmonic coupling coefficient is defined based on the established HSS model. This is used to reveal the harmonic coupling mechanism of the multi-converter topology, and analyze the harmonic interaction characteristics of the cascade and parallel converters. Then, the harmonic coupling coefficient is used to quantify the influence of key parameters such as filter inductor, current loop and phase-locked loop on the harmonic interaction of MGCCS. Finally, the results of the HSS model, Matlab/Simulink model and RT-LAB model are compared to verify the accuracy of the proposed HSS model and the effectiveness of the harmonic coupling coefficient theory applied to the AC/DC harmonic interaction analysis of the system.This work is supported by the National Natural Science Foundation of China (No. 51977127).Key words: multiple grid-connected-converter system; AC/DC harmonic interaction; harmonic state space;harmonic coupling coefficient; grid-connected current0 引言在“双碳”目标的背景下,以高渗透率的可再基金项目:国家自然科学基金项目资助(51977127);上海市科学技术委员会项目资助(19020500800);上海市教育发展基金会和上海市教育委员会“曙光计划”项目资助(20SG52) 生能源、高比例的电力电子设备、高速增长的新型负荷为主要特征的新型电力系统逐步形成[1-4]。
Q uenched Mull ins 2Herr i ng 生长模型的反常动力学标度陈庆虎1 ,2 , 焕1刘 ( 1 . 浙江大学 物理系 ,浙江 杭州 310027 ;2 浙江师范大学 海峡两岸统计物理与凝聚态理论研究中心 ,浙江 金华 321004)摘 要 :利用中间时刻 r z 0 t 0 L z ( L 是系统尺寸且 r 0 L ) 等时的界面高度关联函数 G ( r , t ) , 通 过直接数值积分 ,我们研究了 Q u enched Mullins 2Herring 生长模型的反常动力学标度行为 ,并得到 了生长界面的局域粗糙指数 αloc = 0 . 97 ( 1) , 接近于 Mullins 2Her ring 生长模型的理论值 1 ,我们还 计算了表征反常标度行为的临界指数 θ = 0 . 51 ( 3) 和局域生长指数 β = 0 . 29 ( 2) . 关键词 :生长模型 ;反常动力学标度 ;异常粗糙 中图分类号 :O469 ;O488文献标识码 :A文章编号 :1001 - 2443 (2009) 04 - 0307 - 04引 言生长界面的动力学问题是凝聚态物理和材料科学的重要研究内容. 界面的生长过程表现出普适的动力 学标度行为 ,这种普适性可用一系列的临界指数来表征. 为了研究不同的生长过程 ,人们建立了很多生长模 型 ,如固体 - 固体模型 ( solid 2o n 2solid m o del ,简称 SO S 模型) . 简化了的 SO S 模型常被广泛地用来研究界面的平衡态和非平衡态的性质 . 著名的 Kardar - Parisi - Zhang ( KPZ ) 方程1 ,2能很好的描述 SO S 模型在某些条 件下非平衡态的生长过程 . 在这个生长过程中 ,每个粒子的扩散或吸收是瞬时发生的 ,并且这个短暂的过程 只允许出现一次. 而在晶体的分子束外延 ( m olecular 2bea m 2epita xy ,简称 MB E ) 生长过程中 ,粒子的扩散遵从 A rr h enius 型运动规律 . 与 SO S 模型不同 ,在 MB E 生长中 ,粒子总是以一定的几率运动 ,它的运动受到两个R 0 e - E ( x ) / k T因素的影响 :一是粒子的跳跃几率 R , 另一个是粒子的沉积几率 R d . 粒子以几率 R = 跳动 , 其B 中激活能 E ( x ) = E 0 + n E b i n d 与位置有关 , E 0 是自由粒子的激活能 , n 是粒子在 x 处的最近邻化学键个数 ,E bi n d 是每个化学键的束缚能 ,而沉积几率则表示为单位时间内每一位置沉积的粒子数. 从上述分析我们可以看出 , KPZ 方程已不能用来描述 MB E 的生长过程 . 后来 ,人们发现四阶线性连续的 Mullins 2Herring 方程 可以用来描述化学键环境下的 MB E 生长3 - 8 . 在外力存在时 , d = 1 + 1 维的 Mullins 2Herring 模型的方程 表示如下9 ,10 :5 h ( x , t ) 4= - K A h + η( x , t ) + F ( 1)5 t其中 h ( x , t ) 是位置 x 和时间 t 的函数 , 表示生长界面的高度 ; 方程右边第一项用来描述界面在化学势梯度作用下的扩散过程 , 扩散系数为 K ; 第二项是噪音 , 它满足高斯分布〈η( x , h ) 〉= 0 和〈η( x , h ) η( x , h ) 〉= 2 Dδ( x - x )δ( h - h ) , F 是外界驱动力 . 在外力作用下 , 界面出现钉扎 - 脱钉相变 . 在文献 10 中 , 我们利用短时动力学方法研究了 quenched Mullins 2Herring 模型的脱钉相变 . 通过引进 Binder cumulant ,我们计算了总粗糙指数 ,生长指数和动力学指数 . 我们的短时动力学结果与长时稳态的动 力学模拟的结果在误差范围内是一致的9 . 总粗糙指数 α = 1 . 48 ( 4) ,表明界面在生长过程中表现出超常的 粗糙现象 ,即这个模型具有反常动力学标度行为 . 那么表征反常行为的临界指数是怎样的呢 ? 这就是本文重点研究的内容 .1 标度方法和计算公式在很多情况下 , 我们所研究的生长过程虽然远离平衡态 , 但是我们发现 , 类似于平衡态的二级相变 , 界收稿日期 :2009 - 03 - 15基金项目 :国家自然科学基金 ( 10574107 和 10774128) ; 浙江省自然科学基金重点项目 ( Z 7080203) . 作者简介 :陈庆虎 ( 1966 - ) , 男 , 安徽桐城人 , 教授 、博士生导师 , 主要从事凝聚态理论的研究.面的涨落具有幂律行为 . 生长界面的粗糙程度用总界面高度涨落的均方根 W 表示 , 其定义为 :[ h ( x , t ) - h ( t ) ]2W ( L , t ) =其中 L 是系统格点数 , h ( t ) 是 t 时刻所有格点 x 处 h ( x , t ) 的平均值 . 在临界力 F c , 根据 Family 2Vicsek 假 设 , W 满足动力学标度关系[ 9 - 15 ] :〈 W ( L , t ) 〉= t βf ( L / ξ( t ) )符号〈 〉表示对噪音的统计平均. 因为标度函数 f ( u ) 满足下列关系式 :u α: u 0 1f ( u ) 0const : u 0 1 再加上关联函数为 ξ( t ) ~ t 1/ z , 重写上述标度关系 :t β: t 0L α : t 0Lz Lz〈 W ( L , t ) 〉~这里α,β = α/ z 和 z 分别是总粗糙指数 , 生长指数和动力学指数 . 我们已在短时临界动力学的研究中获得 :α= 1 . 48 ( 4) ,β = 0 . 835 ( 1) , z = 1 . 77 ( 5)[ 10 ]. 我们知道 , 总粗糙指数 α 表征界面达到饱和时的形貌. 因为〈 W 〉sa t ≡ W ( L , t → ∞) ~ L α , 在热力学极限下 , 当 α < 1 时〈, W 〉sat / L 0 Lα- 1趋于零 , 表明界面具有自相 似性 ; 当 α > 1 〈, W 〉sa t / L 0 Lα- 1发散 ,说明界面经过很长一段时间的生长后 ,形貌变化很大 ,与初始的形貌 不再相同 ,系统表现出超常的界面粗糙现象 .通常的 Family - Vicsek 动力学标度存在的前提是界面生长具有自相似性 ,系统没有特征长度 (除了系统尺寸) ,标度变换不会改变系统的物理性质 . 当界面存在超常粗糙行为时 ,不是所有的长度尺寸都等价 ,即局 域的 (l ocal ) 界面涨落和系统总的界面涨落满足不同的动力学标度. 界面局域涨落的标度特性可以用等时高 度关联函数来描述 ,其定义为11 ,12 ,16 ,17:G ( r , t ) =〈[ h ( x + r , t ) - h ( x , t ) ]2〉-其中 r 0 L , 直线符号 表示对所有格点的平均〈, 〉仍是对噪音的统计平均 . 界面局域涨落的标度形式为 : 〈 G ( r , t ) 〉= r 2αf A ( r/ ξ( t ) )f A ( u ) 满足下列标度关系 :u - 2(α- αlo c): u 0 1f A ( u ) 0u - 2α: u 0 1当系统达到饱和以前 ,ξ( t ) - t 1/ z , 当系统达到饱和以后 ,ξ( t ) ~ L , 即 :关联长度被系统尺寸截断. 所以 , 关联长度变为 ξ( t ) ~ L g ( t 1/ z/ L ) , 而 g ( x ) 满足 :x : x 0 1g ( x ) ~const : x 0 1所以得到高度关联函数的标度关系为 :r2αlo c t 2β: r z0 G ( r , t ) ~ r 2αl oc L 2θ: r z 0 L zt 0 L z 0 tt 2β : t 0 r z0 Lzα- αlo c其中 ,θ = α - αlo c 是一个新的指数 , 表征界面的反常生长行为 ,β = θ=是局域生长指数. 由此可zz见 , 通过测量在某中间时刻 ( r z 0 t 0 L z) 的关联函数 , 我们就可以得到界面的局域粗糙指数 αlo c ;根据已知的 总粗糙指数 α和动力学指数 z [ 10 ] , 就能分别得到 θ,β 的值.在目前的模拟中 , 我们取系统尺寸 L = 16384 , 外界驱动力在临界区域 F = 0 . 9877 [ 10 ] . 与以前的模拟过程类似[ 9 , 10 ], 我们假定公式 ( 1) 中系数 K = 1 , 时间间隔 Δt = 0 . 01 ( 任意时间单位) 和噪音中 D = 0 . 5 , 高度h ( x , t ) 为任意长度单位 , 利用下面的迭代公式更新每个格点处界面的高度 :∧h ( x , t + Δt ) = h ( x , t ) + Δt { F + η( x , h ) -- 4 h ( x + 1 , t ) + h ( x + 2 , t ) }h ( x - 2 , t ) - 4 h ( x - 1 , t ) + 6 h ( x , t )2其中 , 根据 ( 2 D ) 2=η ,η(x , h ) ∈ [ -3 , 3 , h 是 h ( x , t ) 的整数部分 .∧∧3计算结果和分析图 1 画出了在五个不同时刻的生长界面. 从图中 , 我们可以看出 , 在 t = 100 ( 任意时间单位) 时 , 界面比2 较光滑 , 随着时间的推移 , 如 t = 500 , 1000 , 2000 ( 任意时间单位) 时的形貌所示 , 界面变得越来越粗糙 , 表明界面没有自相似性 , 其生长过程具有反常动力学特性 .图 1 在 t = 100 , 500 , 1000 , 2000 ( 由下到上) 个时间单位时 , 图 2 在 t = 500 , 1000 , 2000 ( 由下到上) 个时间单位时 ,体系中 1000 个格点的界面形貌 ( 任意长度单位) 高度的关联函数〈 G ( r , t ) 〉与 r 的双对数关系曲线 , 虚线的斜率为 1 . 94 .图 2 是关联函数〈 G ( r , t ) 〉和 r 在时间 t = 500 , 1000 , 2000 ( 任意时间单位) 时的双对数曲线. 从图中可以看出 , 当 r 比 较 小 时 , 三 条 曲 线 G ( r , t ) 和 r 均 呈 现 出 幂 律 关 系 , 对 应 于 上 述 标 度 关 系 中 的 条 件 : r z 0 t 0 L z, 在同一时刻 , 有 l o g G ( r , t ) ~ 2αlo c lo g r . 我们还观察到 , 在不同时刻 , 三条曲线并未重合 , 说明界 面的局域生长行为是时间的函数 , 与图 1 相一致. 图中三个不同时刻的直线部分的斜率约为 1 . 94 , 如图中虚线所示 , 即 2αloc = 1 . 94 ( 2) , 由此得到 αlo c = 0. 97 ( 1) , 接近于 ,与 Mullins 2Herring 模型中的结果一致11 ,13. 根据θ = α - αlo c ,α = 1 . 48 ( 4) 我们得到了反常临界指数θ = 0 . 51 ( 3) ; 再根据β = θ和 z = 1 . 77 ( 5) , z 得到了局域生长指数 β = 0 . 29 ( 2) .3结 论通过以前计算系统总的界面涨落 , 我们得到了总粗糙指数 α = 1 . 48 ( 4) , 这是系统表现出反常动力学标度的信号 . 利用等时的高度关联函数 G ( r , t ) , 我们计算了中间时刻 r z 0 t 0 L z的 G ( r , t ) ~ r 关系 , 得到了界面局域的粗糙指数 αlo c = 0 . 97 ( 1) , 接近于 1 . 利用标度关系式 θ = α - αlo c以及β = θ, 我们得到了表征 z界面反常生长行为的临界指数 θ = 0 . 51 ( 3) , 同时也得到了局域生长指数 β = 0 . 29 ( 2) .参考文献 :1K ARDAR M , PAR ISI G , ZHAN G Y - C. Dynamic Scaling of G rowing Int erfaces J . Phys Rev L et t , 1986 , 56 : 889 ; Medina E , H w a T , K ardar M and Zhang Y - C. Burgers equatio n wit h c o rrelat ed noise : Reno r malizatio n 2gro up analysis and applicatio ns to direct ed polymers a nd int erface growt h J . Phys Rev A , 1989 ,39 :3053 . KIM J M , K OS T ERL IZ J M . Growt h in a rest rict ed solid 2o n 2solid mo del J . Phys Rev L et t , 1989 ,62 :2289 ; Amar J G and Family F , Phys. Phase t ransitio n in a rest rict ed solid 2o n 2solid surface 2growt h mo del in 2 + 1 di mensio ns J . Rev L et t , 1990 ,64 :543 ; G uo H , G ro ssmann B , a nd G rant M . K inetics of int erface growt h in driven syst ems J . Phys Rev L et t , 1990 ,64 :1262 .R # CZ Z , PL ISCH KE M . Widt h dist ributio n f o r ( 2 + 1) 2di mensio nal growt h and depo sitio n p rocesses J . Phys Rev E , 1994 ,50 :3530 . SARMA S D , TAMBOR EN EA P I. A new universalit y class f o r kinetic growt h : One 2di mensio nal molecular 2beam epit axy J . Phys Rev L et t , 1991 ,66 :325 .L A I Z W , SARMA S D. K inetic growt h wit h surface rel axatio n : Co n tinuum versus ato mistic mo dels J . Phys Rev L et t , 1991 ,66 :2348 . SARMA S D , GHA I SAS S V . St ep bunching as a cha otic pat t ern fo r matio n p r ocess J . Phys Rev L et t , 1992 ,69 :3758 .W I LB Y M R , VV ED E NS KY D D , ZAN G W I LL A. Scaling in a solid 2o n 2solid mo del of epit axial growt h J . Phys Rev B , 1992 ,46 :12896 . KIM J M , SARMA S D. G rowt h in a rest rict ed 2curvat ure mo d el J . Phys Rev E , 1993 ,48 :2599 .23 4 5 6 7 8L EE C , KIM J M . Depinning t ransitio n of t he Mullins 2H erring equatio n wit h an ext ernal driving f o rce and quenched rando m diso r der J . Phys Rev E , 2006 ,73 :016140 .L IU H , ZHOU W , N I E Q - M , CHEN Q - H. Depinning t r ansitio n of t he quenched Mullins 2H erring equatio n :A sho rt 2time dynamic met ho d J . Phys L et t . A , 2008 ,372 :7077 .FORR ES T B M , TAN G L H. Surface ro ughening in a hypercube 2st acking mo d el J . Phys Rev L et t , 1990 ,64 :1405 .L δ P EZ J M , RODR I GU EZ M A , CU ERNO R. Superro ughening versus int r insic ano malo us scaling of surfaces J . Phys Rev E , 1997 ,56 : 3993 .L δ P EZ J M . Scaling app r oach to calculat e critical expo nent s in ano m alo u s surface ro ughening J . Phys Rev L et t , 1999 ,83 :4594 . L EE J H , KIM S K , KIM J M . G rowt h wit h surface curvat u re o n quenched pot entials J . Phys Rev E , 2000 ,62 :3299 .RAMSCO J J , L δ P E Z J M , RODR I GU EZ M A. Int erface depinning in t he absence of an ext ernal driving f o rce J . Phys Rev E , 2001 ,64 : 066109 .D IC K MAN R , MU ΣOZ M A. Int erface scaling in t he c o nt act p r ocess J . Phys RevE , 2000 ,62 :7632 .K ARA BACA K T , ZHAO Y - P , WAN G G - C , L U T - M . G rowt h f r o nt ro ughening in silic o n nit ride f il ms by plasma 2enhanced chemical vapo r de po sitio n J . Phys Rev B , 2002 ,66 :075329 .9 10 11 12 13 14 15 16 17Anomal o us Dyna m ic Scal i ng of Super 2R oughen ing in the Q uenc hedMul l i ns 2 H erring G ro wth ModelCH EN Qing 2hu1 ,2, L IU Huan1(1 . Depart ment of Physics , Zhe jiang U n iversit y , H angzho u 310027 , China ; 2 . Cent er f o r St atistical and The o retical Co ndensed Mat t er Physics ,Zhejiang No r mal U n iversit y , J i nhua 321004 , China )Abstract : By measuring equal 2time height 2height co r relati o n f u ncti o n G ( r , t ) fo r inter m ediate times r z0 t 0 L z( L is t h e system size and r 0 L ) , we st u dy t h e ano m al o u s dynamic scaling of super - ro u ghening in t h equenched Mullins - Her ring growt h m o del in a direct numerical integrati o n scheme . The local r o u ghnessexpo nent αlo c = 0 . 97 ( 1) is evaluated , cl o se to t he unit w hich is o btained t heo retically in t he Mullins 2H erri ng growt h m o del . A new e x po nent θ = 0 . 51 ( 3) , characterizing t he ano maly of t he growt h surf ace , and t he l o calgrowt h e x po n ent β = 0 . 29 ( 2) are also ext racted.K ey w ords :growt h m o d el ; ano m al o u s dynamic scaling ; super 2ro u ghness。
a r X i v :c o n d -m a t /0112338v 1 [c o n d -m a t .s t a t -m e c h ] 18 D e c 2001System size resonance in coupled noisy systems and in the Ising modelA.Pikovsky and A.ZaikinDepartment of Physics,University of Potsdam,Postfach 601553,D-14415Potsdam,GermanyM. A.de la CasaDepartamento F´ısica Fundamental,Universidad Nacional de Educaci´o n a Distancia,28040Madrid,Spain(Dated:February 1,2008)We consider an ensemble of coupled nonlinear noisy oscillators demonstrating in the thermody-namic limit an Ising-type transition.In the ordered phase and for finite ensembles stochastic flips of the mean field are observed with the rate depending on the ensemble size.When a small periodic force acts on the ensemble,the linear response of the system has a maximum at a certain system size,similar to the stochastic resonance phenomenon.We demonstrate this effect of system size resonance for different types of noisy oscillators and for different ensembles –lattices with nearest neighbors coupling and globally coupled populations.The Ising model is also shown to demonstrate the system size resonance.PACS numbers:05.40.Ca,05.45.-a,05.50+qStochastic resonance has attracted much interest re-cently [1].As was demonstrated in [2],a response of a noisy nonlinear system to a periodic forcing can exhibit a resonance-like dependence on the noise intensity.In other words,there exists a “resonant”noise intensity at which the response to a periodic force is maximally ordered.Stochastic resonance has been observed in numerous ex-periments [3].Noteworthy,the order in a noise-driven system can have a maximum at a certain noise level even in the absence of periodic forcing,this phenomenon being called coherence resonance [4].Being first discussed in the context of a simple bistable model,stochastic resonance has been also studied in complex systems consisting of many elementary bistable cells [5],moreover,the resonance may be enhanced due to coupling [6].In this paper we discuss another type of resonance in such systems,namely the system size res-onance ,when the dynamics is maximally ordered at a certain number of interacting subsystems.Contrary to previous reports of array-enhanced stochastic resonance (cf.also [7]),here we fix the noise strength,coupling,and other parameters;only the the size of the ensemble changes.The basic model to be considered below is the ensemble of noise-driven bistable overdamped oscillators,governed by the Langevin equations ˙x i =x i −x 3i +ε2Dξi (t )+f (t ).(1)Here ξi (t )is a Gaussian white noise with zero mean: ξi (t )ξj (t ′) =δij δ(t −t ′);εis the coupling constant;N is the number of elements in the ensemble,and f (t )is a periodic force to be specified later.In the absence of periodic force the model (1)has been extensively studied in the thermodynamic limit N →∞.It demonstrates an Ising-type phase transition at ε=εc from the disor-dered state with vanishing mean field X =N −1 i x ito the “ferromagnetic”state with a nonzero mean fieldX =±X 0[8].While in the thermodynamic limit the full description of the dynamics is possible,for finite system sizes we have mainly a qualitative picture:in the ordered phase the mean field X switches between the values ±X 0and its average vanishes for all couplings.The rate of switchings depends on the system size and tends to zero as N →∞.For us,of the main importance is the fact that qualita-tively the behavior of the mean field can be represented as the noise-induced dynamics in a potential with one minimum in the disordered phase (at X =0)and two symmetric minima (at X =±X 0)in the ordered phase.Now applying the ideas of the stochastic resonance,one can expect in the bistable case (i.e.in the ordered phase for small enough noise or for large enough coupling)a resonant-like behavior of the response to a periodic ex-ternal force when the intensity of the effective noise is changed.Because this intensity is inverse proportional to N ,we obtain the resonance-like curve of the response in dependence of the system size.The main idea be-hind the system size resonance is that in finite ensembles of noise-driven or chaotic systems the dynamics of the mean field can be represented as driven by the effective noise whose variance is inverse proportional to the sys-tem size [9].This idea has been applied to description of a transition to collective behavior in [10].In [11]it was demonstrated that the finite-size fluctuations can cause a transition that disappears in the thermodynamic limit.Before proceeding to a quantitative analytic descrip-tion of the phenomenon,we illustrate it with direct nu-merical simulations of the model (1),with a forcing term f (t )=A cos(Ωt ).Figure 1shows the linear response function,i.e.the ratio of the spectral component in the mean field at frequency Ωand the amplitude of forcing A ,in the limit A →0.For a given frequency Ωthe depen-dence on the system size is a bell-shaped curve,with a2pronounced maximum.The dynamics of the mean field X (t )is illustrated in Fig.2,for three different system sizes.The resonant dynamics (Fig.2b)demonstrates a typical for stochastic resonance synchrony between the driving periodic force and the switchings of the field be-tween the two stable positions.0.00010.0010.010.1Ω20406080N50010001500linear responseFIG.1:Linear response of the ensemble (1)(D =0.5,ε=2)in dependence on the frequency and the system size N .timemean field X −101a)b)c)FIG.2:The time dependence of the mean field in the en-semble (1)for D =0.5,ε=2,A =0.02,Ω=π/300,and different sizes of the ensemble:(a)N =80,(b)N =35,and (c)N =15.We also depict the periodic force (its amplitude is not in scale)to demonstrate the synchrony of the switchings with the forcing in (b).To describe the system size resonance analytically,we use,following [8],the Gaussian approximation.In this approximation one writes x i =X +δi and assumes that δi are independent Gaussian random variables with zero mean and the variance M .Assuming furthermorethat N −1i δ2i =M and neglecting the odd moments N −1 i δi ,N−1 i δ3i ,as well as the correlations be-tween δi and δj ,we obtain from (1)the equations for X and M :˙X=X −X 3−3MX +Nη(t )+f (t ),(2)1(2+ε)2−24D )and the unstable pointX =0,M =(1−ε+2D(ε−1)2+12D ,b =−0.5+1.5(ε−1)((ε−1)2+12D )−1/2.A better approx-imation valid also beyond a vicinity of the critical pointcan be constructed if we use ¯b =aX −20instead of b ,where the fixed point X 0is given by (4).Having written the en-semble dynamics as a standard noise-driven double-well system (5)(cf.[1,12]),we can use the analytic formula for the linear response R derived in [12].It reads R =NX 20s )s )21+π2Ω2320406080N20406080l i n e a r r e s p o n s eFIG.3:Comparison of the system-size dependencies of the linear response function for frequencies Ω=0.05(circles)and Ω=0.1(squares)with theory (6).The parameters are D =1and ε−εc =2.5(where the the exact εc and the approximate εc =3D are used for the ensemble and the Gaussian approx-imation,respectively).Inlet:Dependence of the system size yielding maximal linear response on the driving frequency Ω(circles:simulations of the ensemble (1),line is obtained by maximizing the expression (6)).It is instructive to compare the response of the noise-driven system (1)with the noise-free case D =0.With-out external force,the ensemble relaxes eventually to asteady state solution with some mean field X ;in this state each oscillator can be in one of stable steady po-sitions of the potential,correspondingly the oscillators form one or two clusters.From the clustering it follows that the linear response does not depend on the number of elements in the ensemble.Our numerical experiments demonstrated also that the response is system-size in-dependent for large forcing amplitudes as well,where,e.g.,the force-induced cluster-mergings occur.Thus,the effect of system size resonance essentially relies on the presence of noise,which breaks the clustering.Above we have considered the system of globally cou-pled nonlinear oscillators (1).The same effect of system size resonance can be observed in a lattice with nearest neighbors coupling as well.In the thermodynamic limit,the Ising-type phase transition occurs in the lattice (if its dimension is larger than one).Similar to the globally cou-pled ensemble,in finite lattices in the ordered phase the switchings between the two stable states of the mean field are observed.With the same argumentation as above we can conclude that the response of the mean field to a pe-riodic forcing can have a maximum at a certain lattice size,while all other parameters (noise intensity,coupling strength,etc.)are kept constant.We illustrate this in Fig.4.As the next example we consider thetwo-dimensionalNr e s p o n s eFIG.4:Filled circles:Response of a two-dimensional lattice of N with nearest-neighbors coupling for A =0.02,T =500,D =0.5,and ε=4.Squares:Response of system (8)(a two-dimensional lattice with D =1.25,ε=30,A =0.1and T =140).Circles:the same as squares,but for a globally coupled lattice with D =1,ε=20,A =0.1and T =100.nearest neighbor Ising model in the presence of a time-dependent external field.The Hamiltonian of the system readsH =−Jijs i s j −A cos(Ωt ) is i ,(7)where J >0and s i =±1.We are interested in the dependence of the response of the mean magnetization m (t )=1Kj(x j −x i )+√40.00010.0010.01Ω20406080100Nlinear responseFIG.5:Linear response (in arbitrary units)of the Ising model (7)for the temperature T =2J slightly below the critical temperature T c =2.269J .above,such transitions occur even in the absence of the additive noise if the system is finite.Thus,the system size resonance should be observed in the lattice (8)as well.We confirm this in Fig.4.Another possible field of application of the system size resonance is the neuronal dynamics (see,e.g.,[18]).Individual neurons have been demonstrated to exhibit stochastic resonance [3,19].While in experiments one can easily adjust noise to achieve the maximal sensitivity to an external signal,it may be not obvious how this ad-justment takes place in nature.The above consideration shows,that changing the number of elements in a small ensemble of coupled bistable elements to the optimum can significantly improve the sensitivity (cf.[5]).More-over,changing its connectivity and/or coupling strength,a neuronal system can tune itself to signals with different frequencies.Concluding,we have shown that in populations of coupled noise-driven elements,exhibiting in the thermo-dynamic limit the Ising-type transition,in the ordered phase (i.e.for relatively small noise and large coupling)the response to a periodic force achieves maximum at a certain size of the system.We demonstrated this effect for the Ising model,as well as for lattices and globally coupled ensembles noisy oscillators.We expect the sys-tem size resonance to occur also in purely deterministic systems demonstrating the Ising-type transition,e.g.in the Miller-Huse coupled map lattice [20].The system size resonance is described theoretically by reducing the dynamics of the mean field to a low-dimensional bistable model with an effective noise that is inverse proportional to the system size.The stochastic resonance in the mean field dynamics then manifests itself as the system size resonance.We thank N.Brilliantov,A.Neiman,J.Parrondo,F.J.de la Rubia and R.Toral for useful discussions and M.Rosenblum for help in analysis of eq.(6).M.C.thanksthe DGESEIC Project No.PB97-0076and the DGES Grantship No.AP9807249358. A.Z.acknowledges fi-nancial support from ESA.[1]L.Gammaitoni,P.H¨a nggi,P.Jung,and F.Marchesoni,Rev.Mod.Phys.70,223(1998);P.Jung,Phys.Reports 234,175(1993).[2]R.Benzi,A.Sutera,and A.Vulpiani,J.Phys.A:Math.,Gen.14,L453(1981).[3]A.Longtin,A.Bulsara,and F.Moss,Phys.Rev.Lett.67,656(1991);A.Simon and A.Libchaber,Phys.Rev.Lett.68,3375(1992);M.L.Spano,M.Wun-Folge,and W.L.Ditto,Phys.Rev.A 46,5253(1992);S.Barbay,G.Giacomelli,and F.Marin,Phys.Rev.E.61,157(2000).[4]A.Pikovsky and J.Kurths,Phys.Rev.Lett.78,775(1997);A.Neiman,P.I.Saparin,and L.Stone,Phys.Rev.E 56,270(1997).[5]P.Jung,U.Behn,E.Pantazelou,and F.Moss,Phys.Rev.A 46,R1709(1992);M.Morillo,J.G´o mez-Ordo˜n ez,and J.M.Casado,Phys.Rev.E 52,316(1995);J.M.Casado and M.Morillo,Phys.Rev.E 52,2088(1995);H.Gang,H.Haken,and X.Fagen,Phys.Rev.Lett.77,1925(1996);P.Jung and G.Mayer-Kress,Phys.Rev.Lett.74,2130(1995).[6]J.F.Lindner et al.,Phys.Rev.Lett.75,3(1995);Phys.Rev.E 53,2081(1996).[7]A.Neiman,L.Schimansky-Geier,and F.Moss,Phys.Rev.E 56,R9(1997);B.Hu and C.Zhou,Phys.Rev.E 61,R1001(2000).[8]R.C.Desai and R.Zwanzig,J.Stat.Phys.19,1(1978).[9]D.Dawson and J.G¨a rtner,Stochastics 20,247(1987);A.S.Pikovsky and J.Kurths,Physica D 76,411(1994);A.Hamm,Physica D 142,41(2000).[10]A.Pikovsky and S.Ruffo,Phys.Rev.E 59,1633(1999).[11]A.S.Pikovsky,K.Rateitschak,and J.Kurths,Z.PhysikB 95,541(1994).[12]P.Jung and P.H¨a nggi,Phys.Rev.A 44,8032(1991).[13]Z.N´e da,Phys.Rev.E 51,5315(1995);J.Javier Breyand A.Prados,Phys.Lett.A 216,240(1996);U.Siewert and L.Schimansky-Geier,Phys.Rev.E 58,2843(1998).[14]M.E.J.Newmann and G.T.Barkema,Monte CarloMethods in Statistical Physics (Clarendon Press,Oxford,1999).[15]C.V.der Broeck,J.M.R.Parrondo,and R.Toral,Phys.Rev.Lett.73,3395(1994);C.V.der Broeck et al.,Phys.Rev.E 55,4084(1997).[16]nda,A.Zaikin,and L.Schimansky-Geier,Chaos,Solitons and Fractals 9,1367(1998).[17]A.Zaikin,J.Kurths,and L.Schimansky-Geier,Phys.Rev.Lett.85,227(2000);A.Zaikin,K.Murali,and J.Kurths,Phys.Rev.E 63,020103(R)(2001).[18]P. A.Tass,Phase Resetting in Medicine and Biol-ogy.Stochastic Modelling and Data Analysis.(Springer-Verlag,Berlin,1999).[19]J.K.Douglas,L.Wilkens,E.Pantazelou,and F.Moss,Nature 365,337(1993);D.F.Russell,L.A.Wilkens,and F.Moss,Nature 402,291(1999).[20]ler and D.A.Huse,Phys.Rev.E 48,2528(1993).。