CN103646152A - 一种基于矩阵指数的电力系统电磁暂态仿真方法 - Google Patents
一种基于矩阵指数的电力系统电磁暂态仿真方法 Download PDFInfo
- Publication number
- CN103646152A CN103646152A CN201310720115.3A CN201310720115A CN103646152A CN 103646152 A CN103646152 A CN 103646152A CN 201310720115 A CN201310720115 A CN 201310720115A CN 103646152 A CN103646152 A CN 103646152A
- Authority
- CN
- China
- Prior art keywords
- matrix
- simulation
- emulation
- state
- vector
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 239000011159 matrix material Substances 0.000 title claims abstract description 126
- 238000000034 method Methods 0.000 title claims abstract description 75
- 238000004088 simulation Methods 0.000 title claims abstract description 69
- 230000001052 transient effect Effects 0.000 title claims abstract description 52
- 239000013598 vector Substances 0.000 claims abstract description 35
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 28
- 238000004458 analytical method Methods 0.000 claims abstract description 20
- 230000003190 augmentative effect Effects 0.000 claims abstract description 6
- 230000008569 process Effects 0.000 claims description 13
- 230000033001 locomotion Effects 0.000 claims description 6
- 230000003416 augmentation Effects 0.000 claims description 5
- 230000008859 change Effects 0.000 claims description 5
- 230000004069 differentiation Effects 0.000 claims description 4
- 230000005284 excitation Effects 0.000 abstract description 4
- 230000009471 action Effects 0.000 abstract description 2
- 230000010354 integration Effects 0.000 description 13
- 230000006870 function Effects 0.000 description 11
- 230000005540 biological transmission Effects 0.000 description 3
- 238000004364 calculation method Methods 0.000 description 3
- 230000008901 benefit Effects 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 238000005183 dynamical system Methods 0.000 description 2
- 230000005611 electricity Effects 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 230000035800 maturation Effects 0.000 description 2
- 230000010355 oscillation Effects 0.000 description 2
- 230000001360 synchronised effect Effects 0.000 description 2
- 238000012360 testing method Methods 0.000 description 2
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 description 1
- 101100499229 Mus musculus Dhrsx gene Proteins 0.000 description 1
- 239000000654 additive Substances 0.000 description 1
- 230000000996 additive effect Effects 0.000 description 1
- 230000003321 amplification Effects 0.000 description 1
- 238000005094 computer simulation Methods 0.000 description 1
- 238000009795 derivation Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 230000005684 electric field Effects 0.000 description 1
- 238000013178 mathematical model Methods 0.000 description 1
- 230000007935 neutral effect Effects 0.000 description 1
- 238000010606 normalization Methods 0.000 description 1
- 238000003199 nucleic acid amplification method Methods 0.000 description 1
- 230000003534 oscillatory effect Effects 0.000 description 1
- 230000004044 response Effects 0.000 description 1
- 230000017105 transposition Effects 0.000 description 1
Images
Landscapes
- Supply And Distribution Of Alternating Current (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
一种基于矩阵指数的电力系统电磁暂态仿真方法,首先在状态分析框架下建立待研究系统的电磁暂态仿真模型:设定仿真步长Δt及相关仿真参数后,启动仿真程序;在每一个仿真步长内,利用当前时刻激励源g(t)、各阶导数信息、状态矩阵A组成增广的状态矩阵采用Scaling and Squaring算法求解矩阵指数与状态变量向量x(t)与p维单位向量组成的增广状态向量进行矩阵向量乘法,结果作为当前时刻的状态变量,由输出方程得到输出向量y(t),写入输出文件,仿真推进一个步长;并在下一个仿真步长内,首先需判断是否有短路、开路和开关动作。本发明能够实现结构复杂且高度不对称的、模型具有刚性特性的电力系统电磁暂态建模仿真。
Description
技术领域
本发明涉及一种用于电力系统电磁暂态仿真建模方法。特别是涉及一种基于矩阵指数的电力系统电磁暂态仿真方法
背景技术
电力系统电磁暂态仿真主要反映系统中电场与磁场的相互影响产生的电气量的变化过程,可得到几十kHz到工频范围内的三相电压电流瞬时值波形。为了准确获取系统动态过程中的高频特性,电磁暂态的仿真步长通常是微秒级。随着直流输电技术的发展与广泛应用,交直流系统互联给电磁暂态仿真提出了新的挑战。近年来我国电力系统实际运行中出现的一些新的仿真场景,如直流输电送端孤岛运行的电磁暂态仿真,系统规模大,仿真时间长,使用传统电磁仿真算法会消耗较多计算资源和计算时间,需要结合模型特性从算法层面提出针对性的改进。
电力系统电磁暂态仿真本质上可归结为对动力学系统时域响应的求取,它包括系统本身的数学模型和与之相适应的数值算法。
当前,电力系统电磁暂态仿真基本框架可分为两类,包括节点分析法(nodal analysis)和状态变量分析法(state-variable analysis)。基于节点分析框架的电磁暂态仿真方法可概括为先采用某种数值积分方法(通常为梯形积分法)将系统中动态元件的特性方程差分化,得到等效的计算电导与历史项电流源并联形式的诺顿等效电路,此时联立整个电气系统的元件特性方程形成节点电导矩阵,如式(1)所示,对其求解即可得到系统中各节点电压的瞬时值。
Gu=i (1)
式(1)所示的节点电导矩阵为线性方程组,可使用各种成熟的线性稀疏矩阵算法库进行求解。节点分析法广泛应用于EMTP、PSCAD/EMTDC等专业的电力系统电磁暂态仿真程序中,工程上也称基于节点分析框架的电磁暂态仿真工具为EMTP类程序。节点分析法的主要优势体现在程序实现难度和仿真计算效率方面,但由于式(1)的节点电导方程本身已将数值积分方法与系统模型融为一体,导致EMTP类程序在求解算法选择方面缺乏灵活性与开放性,同时式(1)已不能给出系统本身的特征信息。
与节点分析法不同,状态变量分析法属于一般性建模方法(general purpose modeling),不仅适于电路与电力系统仿真,同样也适于其它工程领域的动力学系统的建模与仿真。Matlab/SimPowerSystems软件是状态变量分析框架下暂态仿真程序的典型代表。与节点分析框架相比,状态方程在模型的计算求解方面具有高度的开放性和灵活性,可方便地选择与问题相适应的数值积分方法,同时能够提供关于系统各种特征的丰富信息(如系统的特征值),进而能够从全局角度了解系统的动态特性,为各种快速、准确、高效的仿真算法的开发与测试工作提供了便利条件。
应用状态变量分析的基础是形成式(2)所示标准形式的状态-输出方程,此时系统中的电源作为输入u。
在电力系统仿真领域,式(2)可以由改进节点法Modified Nodal Analysis(MNA)通过KCL、KVL等约束关系以及元件伏安特性进行构造而得到MNA模型,再经过一定的正规化处理(regularization)转化而来。MNA模型是形如式(3)的状态-输出方程。
也可以采用一般支路等方法,如Automated State Model Generator(ASMG)方法直接构造得到。基于这些方法得到的电力系统模型能够很容易的与本发明所采用的状态变量分析框架下的电磁暂态仿真程序进行接口。
在数值算法方面,传统数值积分方法可分为显式和隐式两类,不同积分方法所具有的数值稳定性和数值精度各不相同。一般来说,隐式方法处理仿真模型中刚性特征的能力较强。电力系统由于动态过程时间尺度差异较大,系统模型表现出一定刚性,这使得主流电磁暂态软件EMTP类程序采用隐式方法以保证数值稳定性。从计算开销方面来看,隐式方法在每一时步内需求解线性方程组,极大限制了其仿真在大规模系统的能力。与之相对的,传统显式方法无需迭代,在每一时步内的运算量较小,但其有限的数值稳定域使得仿真步长受到约束,综合来看对刚性系统的仿真性能不佳。对于现代电力系统来说,系统中既存在微秒级的电力电子开关动态过程,又存在同步机组的励磁、调速等秒级的机电动态过程,时间尺度差异极大,系统刚性特征十分显著。充分利用状态方程框架在数值算法选择方面的灵活性,结合电力系统电磁暂态仿真的应用场景与特殊需求,发展合适的数值积分方法,是提高电力系统电磁暂态仿真的计算性能和应用前景的关键。
随着电力电子设备如FACTS、HVDC,以及分布式电源并网系统等在现代电力系统内的广泛应用,电磁暂态和机电暂态之间的界限逐渐模糊。由于电磁暂态仿真采用了系统级的详细模型,能准确反映电力电子开关动作,能与现场运行的实际情况较好吻合,其在实际工程应用中显得尤为重要。在一些特定的仿真场景,如直流输电系统对次同步谐振的影响,或由电力电子装置驱动的大惯性机组的运行情况等,通常要求进行数十秒的电磁暂态仿真模拟。这些场景客观上要求电磁暂态仿真采取尽可能大的仿真步长,以缩短仿真计算时间。电力电子设备在电力系统中的广泛应用,使得电力系统动态过程的多时间尺度特性更为明显,仿真模型刚性更强。
发明内容
本发明所要解决的技术问题是,提供一种适合处理电磁暂态模型刚性特性,保证数值稳定性,同时具有较高精度,并可增大仿真步长的数值积分的基于矩阵指数的电力系统电磁暂态仿真方法。
本发明所采用的技术方案是:一种基于矩阵指数的电力系统电磁暂态仿真方法,包括如下步骤:
1)在状态分析框架下建立待研究电力系统的电磁暂态仿真模型,形式为:
y(t)=Lx(t)
其中,A是状态矩阵,g(t)是激励源,L是输出矩阵,x(t)是状态变量,y(t)是输出变量;
2)启动仿真程序,设定仿真时间Tend,仿真步长Δt,多项式近似系数p,设定当前时刻为仿真起始时刻tn=t0,依照仿真需要设置仿真初值xn=x0,其中的下标n表示tn时刻的相关变量;
4)采用Scaling and Squaring算法求解计算矩阵指数
6)由yn=Lxn得到tn时刻的输出向量并写入输出数据文件,更新当前时刻为下一时刻tn=tn+Δt,仿真向前推进一个步长;
7)比较当前时刻tn与仿真时间Tend,判断是否已经抵达仿真结束时刻,若已经达到,则仿真结束;若未达到,则继续到步骤8);
8)判断当前时刻是否有短路和开路系统故障,或开关动作导致系统模型发生改变的事件,若有,重新建立系统仿真模型,更新状态矩阵A、输出矩阵L;若没有,则跳过更新;
9)返回到步骤3)继续进行下一时刻计算,依此迭代进行,直到仿真结束。
步骤3)所述的利用当前时刻激励源及当前时刻激励源各阶导数信息,组成增广的状态矩阵具体方法如下:
计算当前时刻的激励源g(tn),及当前时刻激励源g(tn)的1至p-1阶导数,如果g(t)的解析形式未知,则用数值微分近似tn时刻g(t)的各阶导数,依照W=[g(p-1)(t0),…,g’’(t0),g’(t0),g(t0)]的格式将这些向量组成矩阵W,按照 的格式更新矩阵,其中矩阵J的格式是固定的 其中I是单位矩阵,n是系统维数,n为大于0的整数。
(1)读入算法参数θ,b1,b2,…b13,其中θ是scaling因子选择阈值,bi是Padé逼近系数;
(2)求解scaling因数s,使得s是满足||A/2s||1<θ的最小正整数,其中||·||1是指矩阵的1范数;
(3)对矩阵A进行scaling运算,即A=A/2s;
(4)求解矩阵eA的[13/13]阶Pade有理近似rA,其过程为:
矩阵运算A2=A2,A4=A2A2,A6=A2A4;
矩阵运算U=A[A6(b13A6+b11A4+b9A2)+b7A6+b5A4+b3A2+b1I];
矩阵运算V=A6(b12A6+b10A4+b8A2)+b6A6+b4A4+b2A2+b0I;
由矩阵线性方程(-U+V)rA=U+V求解rA;
本发明的一种基于矩阵指数的电力系统电磁暂态仿真方法,在状态分析框架下,基于矩阵指数运算实现了常微分方程初值问题的求解,能够实现结构复杂且高度不对称的、模型具有刚性特性的电力系统电磁暂态建模仿真,并在大的仿真步长下保证仿真的数值稳定性和数值精度。本发明的方法不仅适于传统电力系统电磁暂态仿真,更由于其具有的良好的刚性处理能力、大步长仿真能力,特别适用于多时间尺度特征显著的现代电力系统的电磁暂态仿真,所以本发明的方法具有很好的适用性。
附图说明
图1是IEEE123节点基准电网算例结构图;
图2是IEEE123节点基准电网算例特征值分布图;
图3是本发明的一种基于矩阵指数的电力系统电磁暂态仿真方法的流程图;
图4是本发明中求解矩阵指数eA所采用的Scaling and Squaring算法流程图;
图5是IEEE123节点算例中配电变压器上流过的C相电流;
图6是IEEE123节点算例中配电变压器上流过的C相电流在区间[0.0332,0.03355]上的放大图,图中可清楚看出本发明提出的矩阵指数方法与隐式梯形法之间的数值精度比较;
图7是IEEE123节点算例中配电变压器上流过的C相电流在区间[0.01664,0.01666]上的放大图,图中可清楚看出本发明提出的矩阵指数方法与4阶显式龙格库塔法之间的数值精度比较;
图8是IEEE123节点算例中配电变压器上流过的C相电流,采用不同数值积分方法时的绝对误差比较。
具体实施方式
在电力系统电磁暂态仿真中,一次系统的主要设备如输电线路、电力变压器等元件可由其相应的等效电路表示,可以采用MNA或ASMG等自动建模方法(本发明采用ASMG方法)建立形如式(4)的电力系统电磁暂态模型。
其中,向量x(t)包含节点电压和一部分支路电流,向量u(t)包含网络中的独立电压源和独立电流源,矩阵C是电感和电容矩阵,G是电导矩阵,B是独立电源的输入矩阵。在采用MNA时,矩阵C可能是奇异的,需要采用额外的正规化方法进行处理。由此即可得到A=-C-1G,g(t)=C-1Bu(t)。此时系统模型成为式(5)形式。
电磁仿真程序中包含的控制系统,如同步电机的励磁、调速、PSS,通常采用控制框图或传递函数等形式建模。它们可以很容易地转化为形如式(5)的状态方程模型。这样就可以由描述一次系统元件连接关系与参数、控制系统结构和参数的输入数据文件,得到可用于状态变量框架下的电磁暂态仿真系统模型。
当系统发生短路、断线故障等导致系统模型发生改变的事件时,应由MNA或ASMG方法重新形成系统的状态模型。
一般来说,根据经典线性系统理论,形如(5)所示的一般性动态系统从给定初始状态(6)开始的运动轨迹,可以使用矩阵指数算子和卷积运算解析地表示为式(7)形式的Volterra积分方程。
x(t0)=x0 (6)
其中,原微分方程组中线性部分,即矩阵A的作用体现在矩阵指数e(t-t0)A一项中,即通常所称的状态转移矩阵。矩阵A的实质含义是状态量之间的相互耦合影响关系,包含了电力系统模型的主要信息。对这一部分的准确求解确保了数值积分算法的数值精度。激励项的作用则体现在卷积积分中,这一项表示状态量受电压源、电流源的直接影响。出于计算复杂度考虑,通过引入族函数定义对卷积项的求解进行转化。族函数的成员定义为:
电力系统仿真中的电压源、电流源通常满足一定的光滑性条件,可以采用p-1次多项式近似激励源g(t),则式(7)重新组织为至这p+1个族函数的加权求和。本发明中选取多项式为p-1阶截断泰勒展开式,即式(10)。将其带入到式(7),记h=t-t0为第n时步仿真步长,xn和xn+1分别是离散时间点tn和tn+1上真值x(tn)和x(tn+1)的数值近似,则结合族函数定义即可得到式(11)。
此时从数学角度来看,已将常微分方程的初值问题求解转化为矩阵指数函数的求解。通过调用成熟、高效的矩阵指数数值算法,(本发明采用Scaling and Squaring方法),每一时步求解p+1个族函数与向量的乘积,即可对该系统进行电磁暂态仿真计算。
为了减少计算量,本发明进一步设计了如下算法:
令W=[g(p-1)(t0),…,g’’(t0),g’(t0),g(t0)]∈Cn×p,并由分块形式组成矩阵
本发明方法通过矩阵指数算子的引入使积分方法具有A稳定性,克服了系统刚性对数值稳定性的限制;同时,通过对系统的线性动态过程进行精确求解,避免了传统数值积分方法所遇到的局部截断误差问题,使算法具有很高的数值精度,这使得仿真可以采用大步长进行计算。本发明所提出算法的性质可以充分满足当前电磁暂态仿真的新需求,且算法的形式可以与MNA方法、ASMG方法等电力系统自动建模方法方便地进行接口,为进一步进行电磁暂态仿真算法的开发奠定了良好的基础。
如图3所示,本发明的一种基于矩阵指数的电力系统电磁暂态仿真方法,具体包括如下步骤:
1)在状态分析框架下建立待研究电力系统的电磁暂态仿真模型,形式为:
y(t)=Lx(t)
其中,A是状态矩阵,g(t)是激励源,L是输出矩阵,上标T是转置运算,x(t)是状态变量,y(t)是输出变量;
2)启动仿真程序,设定仿真时间Tend,仿真步长Δt,多项式近似系数p,设定当前时刻为仿真起始时刻tn=t0,依照仿真需要设置仿真初值xn=x0,其中的下标n表示tn时刻的相关变量;
3)利用当前时刻激励源及当前时刻激励源各阶导数信息,组成增广的状态矩阵具体是指计算当前时刻的激励源g(tn),及当前时刻激励源g(tn)的1至p-1阶导数,如果g(t)的解析形式未知,则用数值微分近似tn时刻g(t)的各阶导数,依照W=[g(p-1)(t0),…,g’’(t0),g’(t0),g(t0)]的格式将这些向量组成矩阵W,按照 的格式更新矩阵,其中矩阵J的格式是固定的 其中I是单位矩阵,n是系统维数,n为大于0的整数。
(1)读入算法参数θ,b1,b2,…b13,其中θ是scaling因子选择阈值,bi是Padé逼近系数;
(2)求解scaling因数s,使得s是满足||A/2s||1<θ的最小正整数,其中||·||1是指矩阵的1范数;
(3)对矩阵A进行scaling运算,即A=A/2s;
(4)求解矩阵eA的[13/13]阶Pade有理近似rA,其过程为:
矩阵运算A2=A2,A4=A2A2,A6=A2A4;
矩阵运算U=A[A6(b13A6+b11A4+b9A2)+b7A6+b5A4+b3A2+b1I];
矩阵运算V=A6(b12A6+b10A4+b8A2)+b6A6+b4A4+b2A2+b0I;
由矩阵线性方程(-U+V)rA=U+V求解rA;
6)由yn=Lxn得到tn时刻的输出向量并写入输出数据文件,更新当前时刻为下一时刻tn=tn+Δt,仿真向前推进一个步长;
7)比较当前时刻tn与仿真时间Tend,判断是否已经抵达仿真结束时刻,若已经达到,则仿真结束;若未达到,则继续到步骤8);
8)判断当前时刻是否有短路和开路系统故障,或开关动作导致系统模型发生改变的事件,若有,重新建立系统仿真模型,更新状态矩阵A、输出矩阵L;若没有,则跳过更新;
9)返回到步骤3)继续进行下一时刻计算,依此迭代进行,直到仿真结束。
下面结合实例,对本发明的一种基于矩阵指数的电力系统电磁暂态仿真方法做出进一步的说明。
以IEEE123节点算例作为实施例,其系统结构如图1所示。该网络电压等级为4.16kV,在节点150处通过变电站与115kV大电网相连,内部设有单相、三相负荷及线路,是一个复杂不对称的辐射状供电网络。通过特征值分析可以发现,系统状态矩阵中含有大量远离实轴的特征值,如图2所示。这一特征在仿真曲线中即表现为高频振荡的动态过程,其频率最高可达1MHz以上,属于典型刚性的高度振荡(highly oscillatory)系统。对于这一刚性系统,传统数值积分方法进行求解会遇到显著的数值精度和稳定性问题。
第一步,在状态分析框架下建立待研究电力系统的电磁暂态仿真模型,其形式为:
y(t)=Lx(t)
对于本实例,采用ASMG进行系统建模,状态向量x(t)包含支路电流和部分节点电压,激励源g(t)是115kV母线的等效电压源;
第二步,启动仿真程序,设定仿真时间Tend,仿真步长Δt,多项式近似系数p。设定当前时刻为仿真起始时刻tn=t0,依照仿真需要设置仿真初值xn=x0。本实例从零状态启动,设定仿真时间0.06s,仿真步长5us,多项式近似阶数为1,即分段线性近似;
第三步,计算当前时刻的激励源g(tn),及其1至p-1阶导数。如果g(t)的解析形式未知,则用数值微分近似tn时刻g(t)的各阶导数。依照W=[g(p-1)(t0),…,g’’(t0),g’(t0),g(t0)]的格式将这些向量组成矩阵W。按照 的格式更新矩阵,其中矩阵J的格式是固定的 本实例中激励源形式已知,即
1)读入算法参数θ,b1,b2,…b13,其中θ是scaling因子选择阈值,bi是Padé逼近系数;
2)求解scaling因数s,使得s是满足||A/2s||1<θ的最小正整数;
3)对矩阵A进行scaling运算,即A=A/2s;
4)求解矩阵eA的[13/13]阶Pade有理近似rA,其过程为:
矩阵运算A2=A2,A4=A2A2,A6=A2A4;
矩阵运算U=A[A6(b13A6+b11A4+b9A2)+b7A6+b5A4+b3A2+b1I];
矩阵运算V=A6(b12A6+b10A4+b8A2)+b6A6+b4A4+b2A2+b0I;
由矩阵线性方程(-U+V)rA=U+V求解rA;
第五步,状态变量向量与p维单位向量 组成增广状态向量,与第四步得到的矩阵指数进行矩阵向量乘法,得到 取结果向量的前n行赋值给xn;
第六步,由yn=Lxn得到tn时刻的输出向量并写入输出数据文件。更新当前时刻为下一时刻tn=tn+Δt,仿真向前推进一个步长;
第七步,比较当前时刻tn与仿真时间Tend,判断是否已经抵达仿真结束时刻。若已经达到,则仿真结束;若未达到,则继续到第八步;
第八步,判断当前时刻是否有短路和开路等系统故障,或开关动作等导致系统模型发生改变的事件,若有,则重新建立系统仿真模型,更新状态矩阵A、输出矩阵L;若没有,则跳过更新;
第九步,回到第三步继续进行下一时刻计算。依此迭代进行,直到仿真结束。
执行仿真计算的计算机硬件环境为Intel Core2Q84002.66GHz CPU,内存容量2GB;软件环境为Windows7操作系统。
分别将本发明所提出的矩阵指数方法(Matrix Exponential)与典型隐式方法——梯形法(Trap)以及典型显式方法——4阶龙格库塔法(RK4)进行比较,并采用有误差控制的变步长Adams方法作为仿真的准确值。对于矩阵指数方法和隐式梯形法,仿真步长设置为5us;龙格库塔法由于数值稳定性限制,需要采用更小的步长才能抑制误差的数值振荡,仿真步长采用0.1us,是其他仿真方法的1/50。
附图5至图8比较了本发明提出的矩阵指数方法和其它几种数值的仿真结果,其中图5为总体图,图6和图7为放大图,图8为绝对误差曲线。
从放大的图6中可以看出,基于矩阵指数的电磁仿真方法所得仿真结果与小步长的RK4方法、以及作为准确值的变步长Adams法较为接近,而梯形法的仿真结果则与其他方法有显著差异,曲线中存在不合理的尖峰,不能真实反映系统较高频率上的振荡特性。
从进一步放大的图7中可以看出,本发明提出的基于矩阵指数的电磁仿真方法所得仿真结果在时轴的每一个离散点,即图7中的黑色圆圈,相比较于小步长的龙格库塔法,与准确值更靠近。
图8给出了不同仿真方法所得曲线与准确值之间绝对误差的比较,图中的纵坐标是10为底的对数坐标。本发明提出的基于矩阵指数的电磁仿真方法绝对误差在1e-3量级,而传统的隐式梯形法、小步长的显式4阶龙格库塔法绝对误差分别在1e-1、1e-2.5量级。充分验证了本发明提出的仿真方法在数值稳定性和数值精度上的突出优势。
以上算例测试结果证明,本发明提出的一种基于矩阵指数的电力系统电磁暂态仿真方法具有良好的可行性与适用性,为解决电力系统电磁暂态仿真提供了一种很好的解决思路。
Claims (3)
1.一种基于矩阵指数的电力系统电磁暂态仿真方法,其特征在于,包括如下步骤:
1)在状态分析框架下建立待研究电力系统的电磁暂态仿真模型,形式为:
y(t)=Lx(t)
其中,A是状态矩阵,g(t)是激励源,L是输出矩阵,x(t)是状态变量,y(t)是输出变量;
2)启动仿真程序,设定仿真时间Tend,仿真步长Δt,多项式近似系数p,设定当前时刻为仿真起始时刻tn=t0,依照仿真需要设置仿真初值xn=x0,其中的下标n表示tn时刻的相关变量;
5)状态变量向量与p维单位向量 组成增广状态向量,与步骤4)得到的矩阵指数进行矩阵向量乘法,得到 取结果向量的前n行赋值给xn;
6)由yn=Lxn得到tn时刻的输出向量并写入输出数据文件,更新当前时刻为下一时刻tn=tn+Δt,仿真向前推进一个步长;
7)比较当前时刻tn与仿真时间Tend,判断是否已经抵达仿真结束时刻,若已经达到,则仿真结束;若未达到,则继续到步骤8);
8)判断当前时刻是否有短路和开路系统故障,或开关动作导致系统模型发生改变的事件,若有,重新建立系统仿真模型,更新状态矩阵A、输出矩阵L;若没有,则跳过更新;
9)返回到步骤3)继续进行下一时刻计算,依此迭代进行,直到仿真结束。
(1)读入算法参数θ,b1,b2,…b13,其中θ是scaling因子选择阈值,bi是Padé逼近系数;
(2)求解scaling因数s,使得s是满足||A/2s||1<θ的最小正整数,其中||·||1是指矩阵的1范数;
(3)对矩阵A进行scaling运算,即A=A/2s;
(4)求解矩阵eA的[13/13]阶Pade有理近似rA,其过程为:
矩阵运算A2=A2,A4=A2A2,A6=A2A4;
矩阵运算U=A[A6(b13A6+b11A4+b9A2)+b7A6+b5A4+b3A2+b1I];
矩阵运算V=A6(b12A6+b10A4+b8A2)+b6A6+b4A4+b2A2+b0I;
由矩阵线性方程(-U+V)rA=U+V求解rA;
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201310720115.3A CN103646152B (zh) | 2013-12-23 | 2013-12-23 | 一种基于矩阵指数的电力系统电磁暂态仿真方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201310720115.3A CN103646152B (zh) | 2013-12-23 | 2013-12-23 | 一种基于矩阵指数的电力系统电磁暂态仿真方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN103646152A true CN103646152A (zh) | 2014-03-19 |
CN103646152B CN103646152B (zh) | 2016-08-17 |
Family
ID=50251365
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201310720115.3A Active CN103646152B (zh) | 2013-12-23 | 2013-12-23 | 一种基于矩阵指数的电力系统电磁暂态仿真方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN103646152B (zh) |
Cited By (22)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104217074A (zh) * | 2014-08-27 | 2014-12-17 | 天津大学 | 一种基于矩阵指数的电磁暂态隐式降阶仿真方法 |
CN104298809A (zh) * | 2014-08-27 | 2015-01-21 | 天津大学 | 一种基于矩阵指数电磁暂态仿真的非线性建模求解方法 |
CN104376158A (zh) * | 2014-11-05 | 2015-02-25 | 天津大学 | 一种用于矩阵指数的暂态仿真多时间尺度输出方法 |
CN105224754A (zh) * | 2015-10-14 | 2016-01-06 | 清华大学 | 一种基于插值补偿电流开关模型的电力电子仿真方法 |
CN106446428A (zh) * | 2016-09-29 | 2017-02-22 | 全球能源互联网研究院 | 一种开关电路电磁暂态分析方法及分析装置 |
CN107230161A (zh) * | 2017-05-17 | 2017-10-03 | 国网北京市电力公司 | 电力系统仿真算法的评价方法和装置 |
CN107356822A (zh) * | 2017-06-28 | 2017-11-17 | 西安交通大学 | 用于电磁脉冲多端口效应评估的多通道探测系统 |
CN107679287A (zh) * | 2017-09-11 | 2018-02-09 | 三峡大学 | 基于3步4阶隐式泰勒级数法的电磁暂态数值计算方法 |
CN108021999A (zh) * | 2016-11-04 | 2018-05-11 | 中国电力科学研究院 | 一种快速逼近最大负荷功率点的方法及装置 |
CN105468864B (zh) * | 2015-12-14 | 2018-05-15 | 三峡大学 | 基于増维精细积分的高压输电线路电磁暂态数值计算方法 |
CN109086249A (zh) * | 2018-08-02 | 2018-12-25 | 北京知存科技有限公司 | 模拟向量-矩阵乘法运算电路 |
CN110457732A (zh) * | 2019-05-30 | 2019-11-15 | 中国南方电网有限责任公司电网技术研究中心 | 交直流电力系统的混合仿真方法、装置及存储介质 |
CN110738016A (zh) * | 2019-10-12 | 2020-01-31 | 南方电网科学研究院有限责任公司 | 一种电力电子电路暂态仿真插值计算方法 |
CN110968938A (zh) * | 2019-10-31 | 2020-04-07 | 全球能源互联网研究院有限公司 | 一种用于电磁暂态仿真的理想开关过程分析方法及系统 |
CN111697889A (zh) * | 2020-05-06 | 2020-09-22 | 南方电网科学研究院有限责任公司 | 一种基于时域变换的异步电动机仿真建模方法及装置 |
CN112069668A (zh) * | 2020-08-26 | 2020-12-11 | 三峡大学 | 基于微分求积法的电磁暂态快速仿真方法 |
CN112214899A (zh) * | 2020-10-16 | 2021-01-12 | 哈尔滨理工大学 | 双轴励磁同步发电机的2s-dirk电磁暂态建模方法 |
CN112784500A (zh) * | 2021-03-22 | 2021-05-11 | 重庆邮电大学 | 基于深度学习和pscad的电磁暂态仿真模型的敏捷生成方法 |
CN112989739A (zh) * | 2021-04-20 | 2021-06-18 | 北京华大九天科技股份有限公司 | 一种Trap-Gear时间离散格式的时间步长设定方法 |
CN113128072A (zh) * | 2021-05-13 | 2021-07-16 | 清鸾科技(成都)有限公司 | 传递函数高精度仿真方法、装置、存储介质及电子设备 |
CN113420433A (zh) * | 2021-06-18 | 2021-09-21 | 中国科学院电工研究所 | 对等控制方式下低压交直流系统可扩展建模及分析方法 |
CN115238465A (zh) * | 2022-06-24 | 2022-10-25 | 南方电网科学研究院有限责任公司 | 电磁暂态仿真中列降阶模型的执行时间计算方法和装置 |
-
2013
- 2013-12-23 CN CN201310720115.3A patent/CN103646152B/zh active Active
Non-Patent Citations (7)
Title |
---|
ALTAN ODABASIOGLU等: "PRIMA:Passive Reduced-Order Interconnect Macromodeling Algorithm", 《IEEE TRANSACTIONS ON COMPUTER-AIDED DESIGN OF INTEGRATED CIRCUITS AND SYSTEMS》 * |
AWAD H.AL-MOHY等: "Computing the Action of the Matrix Exponential, with an Application to Exponential Integrators", 《SIAM JOURNAL ON SCIENTIFIC COMPUTING》 * |
JEAN MAHSEREDJIAN等: "Simulation Tools for Electromagnetic Transients in Power Systems:Overview and Challenges", 《IEEE TRANSACTIONS ON POWER DELIVERY》 * |
MARLIS HOCHBRUCK等: "Exponential integrators", 《ACTA NUMERICA》 * |
NICHOLAS J.HIGHAM: "The Scaling and Squaring Method for the Matrix Exponential Revisited", 《SIAM JOURNAL ON MATRIX ANALYSIS AND APPLICATIONS》 * |
O.WASYNCZUK等: "Automated State Model Generation Algorithm for Power Circuits and Systems", 《IEEE TRANSACTIONS ON POWER SYSTEMS》 * |
李鹏等: "基于Krylov子空间的大规模配电网络模型整体化简方法", 《电网技术》 * |
Cited By (37)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104217074B (zh) * | 2014-08-27 | 2017-04-12 | 天津大学 | 一种基于矩阵指数的电磁暂态隐式降阶仿真方法 |
CN104298809A (zh) * | 2014-08-27 | 2015-01-21 | 天津大学 | 一种基于矩阵指数电磁暂态仿真的非线性建模求解方法 |
CN104217074A (zh) * | 2014-08-27 | 2014-12-17 | 天津大学 | 一种基于矩阵指数的电磁暂态隐式降阶仿真方法 |
CN104376158B (zh) * | 2014-11-05 | 2017-05-10 | 天津大学 | 一种用于矩阵指数的暂态仿真多时间尺度输出方法 |
CN104376158A (zh) * | 2014-11-05 | 2015-02-25 | 天津大学 | 一种用于矩阵指数的暂态仿真多时间尺度输出方法 |
CN105224754A (zh) * | 2015-10-14 | 2016-01-06 | 清华大学 | 一种基于插值补偿电流开关模型的电力电子仿真方法 |
CN105224754B (zh) * | 2015-10-14 | 2018-08-10 | 清华大学 | 一种基于插值补偿电流开关模型的电力电子仿真方法 |
CN105468864B (zh) * | 2015-12-14 | 2018-05-15 | 三峡大学 | 基于増维精细积分的高压输电线路电磁暂态数值计算方法 |
CN106446428A (zh) * | 2016-09-29 | 2017-02-22 | 全球能源互联网研究院 | 一种开关电路电磁暂态分析方法及分析装置 |
CN108021999A (zh) * | 2016-11-04 | 2018-05-11 | 中国电力科学研究院 | 一种快速逼近最大负荷功率点的方法及装置 |
CN108021999B (zh) * | 2016-11-04 | 2022-02-22 | 中国电力科学研究院 | 一种快速逼近最大负荷功率点的方法及装置 |
CN107230161A (zh) * | 2017-05-17 | 2017-10-03 | 国网北京市电力公司 | 电力系统仿真算法的评价方法和装置 |
CN107356822B (zh) * | 2017-06-28 | 2019-06-25 | 西安交通大学 | 用于电磁脉冲多端口效应评估的多通道探测系统 |
CN107356822A (zh) * | 2017-06-28 | 2017-11-17 | 西安交通大学 | 用于电磁脉冲多端口效应评估的多通道探测系统 |
CN107679287A (zh) * | 2017-09-11 | 2018-02-09 | 三峡大学 | 基于3步4阶隐式泰勒级数法的电磁暂态数值计算方法 |
CN107679287B (zh) * | 2017-09-11 | 2021-07-13 | 三峡大学 | 基于3步4阶隐式泰勒级数法的电磁暂态数值计算方法 |
CN109086249A (zh) * | 2018-08-02 | 2018-12-25 | 北京知存科技有限公司 | 模拟向量-矩阵乘法运算电路 |
CN109086249B (zh) * | 2018-08-02 | 2023-10-20 | 北京知存科技有限公司 | 模拟向量-矩阵乘法运算电路 |
CN110457732A (zh) * | 2019-05-30 | 2019-11-15 | 中国南方电网有限责任公司电网技术研究中心 | 交直流电力系统的混合仿真方法、装置及存储介质 |
CN110457732B (zh) * | 2019-05-30 | 2023-02-28 | 南方电网科学研究院有限责任公司 | 交直流电力系统的混合仿真方法、装置及存储介质 |
CN110738016A (zh) * | 2019-10-12 | 2020-01-31 | 南方电网科学研究院有限责任公司 | 一种电力电子电路暂态仿真插值计算方法 |
CN110968938A (zh) * | 2019-10-31 | 2020-04-07 | 全球能源互联网研究院有限公司 | 一种用于电磁暂态仿真的理想开关过程分析方法及系统 |
CN110968938B (zh) * | 2019-10-31 | 2024-03-15 | 全球能源互联网研究院有限公司 | 一种用于电磁暂态仿真的理想开关过程分析方法及系统 |
CN111697889A (zh) * | 2020-05-06 | 2020-09-22 | 南方电网科学研究院有限责任公司 | 一种基于时域变换的异步电动机仿真建模方法及装置 |
CN112069668A (zh) * | 2020-08-26 | 2020-12-11 | 三峡大学 | 基于微分求积法的电磁暂态快速仿真方法 |
CN112069668B (zh) * | 2020-08-26 | 2023-06-30 | 三峡大学 | 电磁暂态仿真中基于微分求积法和v变换的矩阵计算方法 |
CN112214899A (zh) * | 2020-10-16 | 2021-01-12 | 哈尔滨理工大学 | 双轴励磁同步发电机的2s-dirk电磁暂态建模方法 |
CN112784500B (zh) * | 2021-03-22 | 2022-07-01 | 重庆邮电大学 | 基于深度学习和pscad的电磁暂态仿真模型的敏捷生成方法 |
CN112784500A (zh) * | 2021-03-22 | 2021-05-11 | 重庆邮电大学 | 基于深度学习和pscad的电磁暂态仿真模型的敏捷生成方法 |
CN112989739B (zh) * | 2021-04-20 | 2021-07-30 | 北京华大九天科技股份有限公司 | 一种Trap-Gear时间离散格式的时间步长设定方法 |
CN112989739A (zh) * | 2021-04-20 | 2021-06-18 | 北京华大九天科技股份有限公司 | 一种Trap-Gear时间离散格式的时间步长设定方法 |
CN113128072A (zh) * | 2021-05-13 | 2021-07-16 | 清鸾科技(成都)有限公司 | 传递函数高精度仿真方法、装置、存储介质及电子设备 |
CN113128072B (zh) * | 2021-05-13 | 2024-01-19 | 清鸾科技(成都)有限公司 | 传递函数高精度仿真方法、装置、存储介质及电子设备 |
CN113420433A (zh) * | 2021-06-18 | 2021-09-21 | 中国科学院电工研究所 | 对等控制方式下低压交直流系统可扩展建模及分析方法 |
CN113420433B (zh) * | 2021-06-18 | 2023-07-21 | 中国科学院电工研究所 | 对等控制方式下低压交直流系统可扩展建模及分析方法 |
CN115238465A (zh) * | 2022-06-24 | 2022-10-25 | 南方电网科学研究院有限责任公司 | 电磁暂态仿真中列降阶模型的执行时间计算方法和装置 |
CN115238465B (zh) * | 2022-06-24 | 2023-04-28 | 南方电网科学研究院有限责任公司 | 电磁暂态仿真中列降阶模型的执行时间计算方法和装置 |
Also Published As
Publication number | Publication date |
---|---|
CN103646152B (zh) | 2016-08-17 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN103646152A (zh) | 一种基于矩阵指数的电力系统电磁暂态仿真方法 | |
CN104298809A (zh) | 一种基于矩阵指数电磁暂态仿真的非线性建模求解方法 | |
CN102184297B (zh) | 适于微网暂态并行仿真的电气与控制系统解耦预测方法 | |
CN104217074B (zh) | 一种基于矩阵指数的电磁暂态隐式降阶仿真方法 | |
CN103700036B (zh) | 一种适于电力系统多时间尺度的暂态稳定性投影积分方法 | |
CN103077268B (zh) | 面向电力系统电磁暂态仿真的状态空间自动建模方法 | |
CN103049617B (zh) | 保留无源性的大规模配电网络电磁暂态仿真模型化简方法 | |
CN103810646B (zh) | 一种基于改进投影积分算法的有源配电系统动态仿真方法 | |
CN103942372A (zh) | 基于fpga的有源配电网暂态实时仿真多速率接口方法 | |
CN103001214B (zh) | 一种基于中性点偏移的配电网三相潮流计算方法 | |
CN105182245A (zh) | 基于无迹卡尔曼滤波的大容量电池系统荷电状态估计方法 | |
CN105574809A (zh) | 基于矩阵指数的电磁暂态仿真图形处理器并行计算方法 | |
CN103440374A (zh) | 基于状态变量节点分析混合框架的电磁暂态仿真建模方法 | |
CN104091092A (zh) | 大规模电力系统小干扰稳定性的特征值分析系统 | |
CN104917174A (zh) | 一种针对主动配电网的静态电压稳定性判定方法 | |
CN105224754A (zh) | 一种基于插值补偿电流开关模型的电力电子仿真方法 | |
CN105353315A (zh) | 一种基于无迹卡尔曼滤波的电池系统荷电状态估计方法 | |
CN101794998A (zh) | 基于多机电力系统单一发电机电磁功率简洁表达方式的暂态稳定在线分析方法 | |
US11569682B2 (en) | System and method for a fast power network simulator | |
CN105182246A (zh) | 基于无迹卡尔曼滤波的并联型电池系统荷电状态估计方法 | |
CN105353314A (zh) | 一种并联型电池系统荷电状态估计方法 | |
CN108054757A (zh) | 一种内嵌无功和电压的n-1闭环安全校核方法 | |
CN104732459A (zh) | 大规模电力系统病态潮流分析系统 | |
CN103678798A (zh) | 一种用于含分布式电源配电网的电磁暂态实时仿真方法 | |
Liu et al. | Fast power system dynamic simulation using continued fractions |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
C14 | Grant of patent or utility model | ||
GR01 | Patent grant |