CN103455668A - 面向状态变量节点分析混合框架的电磁暂态仿真插值方法 - Google Patents

面向状态变量节点分析混合框架的电磁暂态仿真插值方法 Download PDF

Info

Publication number
CN103455668A
CN103455668A CN2013103742583A CN201310374258A CN103455668A CN 103455668 A CN103455668 A CN 103455668A CN 2013103742583 A CN2013103742583 A CN 2013103742583A CN 201310374258 A CN201310374258 A CN 201310374258A CN 103455668 A CN103455668 A CN 103455668A
Authority
CN
China
Prior art keywords
delta
equation
constantly
combination frame
state
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
Application number
CN2013103742583A
Other languages
English (en)
Other versions
CN103455668B (zh
Inventor
洪潮
周保荣
姚文峰
王成山
李鹏
于浩
宋关羽
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Tianjin University
Research Institute of Southern Power Grid Co Ltd
Original Assignee
Tianjin University
Power Grid Technology Research Center of China Southern Power Grid Co Ltd
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Tianjin University, Power Grid Technology Research Center of China Southern Power Grid Co Ltd filed Critical Tianjin University
Priority to CN201310374258.3A priority Critical patent/CN103455668B/zh
Publication of CN103455668A publication Critical patent/CN103455668A/zh
Application granted granted Critical
Publication of CN103455668B publication Critical patent/CN103455668B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

一种面向状态变量节点分析混合框架的电磁暂态仿真插值方法：基于节点方程和状态方程分别建立研究系统的电磁暂态仿真模型；对由节点方程和状态方程构成的混合框架进行由t时刻到t+Δt时刻的梯形法积分求解；检测所有开关动作情况；改变开关状态，重新形成暂态计算矩阵；重新检测所有开关动作情况；运用线性插值方法，对混合框架中的节点方程和状态方程ts和时刻值外插到
Figure DDA0000371483660000012
时刻值；对混合框架进行由
Figure DDA0000371483660000013
时刻到
Figure DDA0000371483660000014
时刻采用两步半步长后向欧拉法进行重新初始化；比较时刻与时刻t+Δt的前后关系，若
Figure DDA0000371483660000016
则进行一步线性插值得到t+Δt时刻的值，否则对混合框架进行由
Figure DDA0000371483660000017
时刻到
Figure DDA0000371483660000018
时刻的梯形法积分求解。本发明实现了开关动作时刻系统状态的计算与重初始化，消除了开关动作之后数值振荡问题。

Description

面向状态变量节点分析混合框架的电磁暂态仿真插值方法
技术领域
本发明涉及一种电力系统电磁暂态仿真的方法。特别是涉及一种面向状态变量节点分析混合框架的电磁暂态仿真插值方法。
背景技术
电力系统电磁暂态仿真本质上可归结为对动力学系统时域响应的求取，它包括系统本身的数学模型和与之相适应的数值算法，其基本仿真框架可分为两类，即节点分析法（NodalAnalysis）和状态变量分析法（State-space Analysis）。基于节点分析的电磁暂态仿真方法可概括为先采用某种数值积分方法（通常为梯形积分法）将系统中动态元件的特性方程差分化，得到等效计算电导与历史项电流源并联形式的诺顿等效电路，此时联立整个电气系统的元件特性方程形成节点电导矩阵，如式（1）所示，求解得到系统中各节点电压的瞬时值。
Gu=i   （1）
式（1）所示的节点电导矩阵为线性方程组，可使用各种成熟的线性稀疏矩阵算法库进行求解。节点分析法广泛应用于EMTP、PSCAD/EMTDC等专业的电力系统电磁暂态仿真程序中，工程上也称基于节点分析框架的电磁暂态仿真工具为EMTP类程序。节点分析法的主要优势体现在程序实现难度和仿真计算效率方面，但由于式（1）所示的节点方程本身已将数值积分方法与系统模型融为一体，导致EMTP类程序在求解算法选择方面缺乏灵活性与开放性，同时式（1）亦不能给出系统本身的特征信息。
与节点分析法不同，状态变量分析法属于一般性建模方法（general purpose modeling），不仅适于电路与电力系统仿真，同样也适于其它形式的动力学系统的建模与仿真。Matlab/SimPowerSystems软件是状态变量分析框架下暂态仿真程序的典型代表。应用状态变量分析的基础是形成式（2）所示标准形式的状态输出方程。上述标准形式的状态方程可通过状态空间自动建模方法ASMG（Automated State Model Generation）得到，并使用各种显式或隐式积分方法进行求解。
x · = Ax + Bu y = Cx + Du - - - ( 2 )
与节点分析法相比，状态方程在模型的计算求解方面具有高度的开放性和灵活性，同时能够提供关于系统特征的丰富信息（如系统的特征值），进而能够从全局角度了解系统的动态特性，为各种快速、准确、高效的仿真算法的开发与测试工作提供了便利条件。但其缺点同样非常明显，与节点法相比，状态变量分析法的建模过程要复杂的多，同时在处理开关动作方面效率较低。因此，基于节点分析框架实现电磁暂态仿真程序成为当前的主流方法。然而，随着近年来数学与系统科学研究的不断发展，出现了为数众多的新方法、新工具值得为电力系统仿真研究所借鉴。状态变量分析作为一种一般性建模方法，与其他领域研究成果的联系更为紧密，在研究方法的衔接与共享方面具有节点法难以比拟的优势。因此，有研究学者在传统节点法仿真框架的基础上，进一步将线性非时变状态方程的求解融入其中，提出了状态变量节点分析混合仿真框架（Combined State Space Nodal Method），并实现了两部分的联立求解，从而在一定程度上解决了两种仿真框架的兼容问题，扩展了传统电磁暂态仿真程序的建模、仿真与分析能力。
在状态变量节点分析混合仿真框架中，节点方程与状态方程的接口可以选择联系节点电压和联系节点注入电流，如图1所示。其中，如图2所示的电感电阻串联电路，整理成状态方程的形式时是以联系节点电压e为输入量，联系节点注入电流i为输出量；而如图3所示的电容电阻并联电路，整理成状态方程的形式时是以联系节点注入电流i为输入量，联系节点电压e为输出量。考虑电力系统实际，本发明选取联系节点电压e为输入量，联系节点注入电流i为输出量。
在混合仿真过程中，节点方程（Gu=i）的形成与处理方式与传统节点分析法相同，而状态方程部分 x · = Ax + Bu y = Cx + Du 则首先建立其梯形法差分化形式：
x ( t + Δt ) = A ^ x ( t ) + B ^ e N ( t ) + B ^ e N ( t + Δt ) - - - ( 3 )
iN(t+Δt)=Cx(t+Δt)+DeN(t+Δt)   （4）
其中Δt为仿真步长，并有 A ^ = ( I - Δt 2 A ) - 1 ( I + Δt 2 A ) , B ^ = ( I - Δt 2 A ) - 1 ( Δt 2 B ) , I为单位矩阵。
联立式（3）与式（4）将状态变量x(t+Δt)消去后便可得到状态方程的输入输出关系：
iN(t+Δt)=iN_hist+WseN(t+Δt)   （5）
其中，iN_hist是由上一步的值计算得到，并有 i N _ hist = C ( A ^ x ( t ) + B ^ u ( t ) ) , W s = C B ^ + D . 式（5）可视为等效节点电导矩阵Ws和等效节点注入电流源列向量iN_hist的并联，具有与节点方程中多相耦合元件差分后的特性方程相似的形式，如图4所示。由于状态方程的联系节点必然包含在节点方程中，因此可将式（5）中的等效节点电导矩阵Ws和等效节点注入电流源列向量iN_hist按照联系节点编号分别插入到已有的节点方程中，从而得到整个系统的计算方程：
Gglobalu=iglobal   （6）
利用式（6）即可计算系统各节点电压的瞬时值，从而实现了状态方程与节点方程的混合求解。
由于上述状态变量节点分析混合框架是在传统节点分析框架基础上进行扩展得到的，因此同样存在传统EMTP类程序的一些问题，例如在高压直流输电（HVDC）、柔性交流输电（FACTS）、含分布式电源及储能装置的智能配电网、局部孤岛运行电网的暂态仿真等研究工作中经常需要对各种形式的电力电子装置进行建模求解，此时传统的EMTP类程序受节点分析框架的限制将面临以下挑战：
首先，定步长仿真算法只能在步长的整数倍时刻改变开关状态，这种开关动作时间上的延迟会导致电压电流波形出现不真实的“尖峰”，即非特征谐波。因此，为了保证仿真的准确性，必须要求定步长算法能够精确考虑开关的动作时刻。一种方法是，在出现开关动作的步长内改用小步长积分到开关动作的准确时刻，这样需要重新计算各元件的等效电导且花费更多的时间用于形成节点电导矩阵并求解。一种更为有效且被广泛使用的方式是采用线性插值的方法，此时不用重新进行积分就能“还原”到开关动作时刻开关动作前系统中的各变量值。
其次，除了开关动作延迟引起的非特征谐波问题，在开关动作后使用梯形法还会引起数值振荡问题。
再次，电力电子装置中存在的大量开关模型会在仿真的某一时刻出现由于一个开关动作而引起的连锁反应情况，即同步开关动作（simultaneous switching）。
最后，对于定步长算法，由于电力电子装置的频繁动作，在一个仿真步长内可能会出现多个开关动作时刻，这被称为多重开关（multiple switching）。
因此，在应用状态变量节点分析混合框架实现电力系统电磁暂态仿真、含分布式电源及储能装置的智能配电网、局部孤岛运行电网的暂态仿真时必须要解决状态变量节点分析混合框架处理电力电子开关电路的能力。针对上述问题，传统的基于节点分析的EMTP类程序提出了一系列基于插值算法的处理方法以解决定步长引起的开关动作与仿真步长不一致及多重开关问题，解决同步开关的动作问题，解决开关动作后梯形法引起的数值振荡问题。面向节点分析的插值算法只需对底层的单个元件进行相应方式的插值处理并修正节点方程，即可求出整个系统开关动作后的准确解。而对于状态方程来说，由于其所描述系统的元件组成以及系统结构都是不明确的，无法通过元件级插值进行修正；同时根据所表述系统的不同，其状态变量的组成可能更为丰富，而不仅仅局限于电感电流与电容电压。
基于状态变量节点分析混合框架的电磁暂态仿真方法兼顾了节点分析法与状态变量分析法的优势，提升了传统的基于节点分析框架仿真程序在建模、仿真与分析方面的能力，但由于状态变量节点分析混合框架只是从系统层面解决了节点分析与状态变量分析两种仿真方法的兼容性问题，而针对节点方程的插值方法并不完全适用于混合框架下的状态方程处理，特别是对于混合框架下的状态方程没有考虑到含大规模电力电子电路系统仿真时的插值算法应用问题。因此，需要开发一种面向状态变量节点分析混合框架的电磁暂态仿真插值方法。
发明内容
本发明所要解决的技术问题是，提供一种实现了开关动作时刻开关动作后系统状态的计算与重初始化，消除了开关动作之后数值振荡问题的面向状态变量节点分析混合框架的电磁暂态仿真插值方法。
本发明所采用的技术方案是：一种面向状态变量节点分析混合框架的电磁暂态仿真插值方法，包括如下步骤：
1）基于节点方程和状态方程分别建立研究系统的电磁暂态仿真模型，确定节点方程与状态方程的接口变量，以联系节点电压eN为状态方程的输入量和以状态方程注入联系节点电流iN为状态方程的输出量，即u=eN，y=iN；
2）对由节点方程和状态方程构成的混合框架进行由t时刻到t+Δt时刻的梯形法积分求解；
3）检测所有开关动作情况，如果不存在开关动作则转步骤2），否则通过线性插值的方法，由t和t+Δt时刻节点方程与状态方程中各值内插到最先动作开关ts时刻值，即计算
e N ( t s ) = e N ( t ) + α ( e N ( t + Δt ) - e N ( t ) ) x ( t s ) = x ( t ) + α ( x ( t + Δt ) - x ( t ) ) i N ( t s ) = i N ( t ) + α ( i N ( t + Δt ) - i N ( t ) ) ,
其中插值系数由开关电流最先反向的开关在t和t+Δt时刻电流值is(t)和is(t+Δt)确定，即 α = i s ( t ) i s ( t ) - i s ( t + Δt ) ;
4）改变开关状态，重新形成暂态计算矩阵，对混合框架进行由ts时刻到
Figure BDA0000371483640000043
时刻的半步长后向欧拉法积分求解；
5）重新检测所有开关动作情况，如果没有其他开关动作，则转入下一步，否则返回步骤4）；
6）运用线性插值方法，对混合框架中的节点方程和状态方程ts和
Figure BDA0000371483640000044
时刻值外插到
Figure BDA0000371483640000045
时刻值，此时α=-1.0，即计算 e N ( t s - Δt 2 ) = e N ( t s ) + α ( e N ( t s + Δt 2 ) - e N ( t s ) ) x ( t s - Δt 2 ) = x ( t s ) + α ( x ( t s + Δt 2 ) - x ( t s ) ) i N ( t s - Δt 2 ) = i N ( t s ) + α ( i N ( t s + Δt 2 ) - i N ( t s ) ) ;
7）对混合框架进行由
Figure BDA0000371483640000047
时刻到
Figure BDA0000371483640000048
时刻采用两步半步长后向欧拉法进行重新初始化；
8）比较时刻
Figure BDA0000371483640000049
与时刻t+Δt的前后关系，若
Figure BDA00003714836400000410
则采用如下方法一，否则采用如下方法二；
方法一：进行一步线性插值得到t+Δt时刻的值，此时
Figure BDA00003714836400000411
即计算 e N ( t + Δt ) = e N ( t s ) + α ( e N ( t s + Δt 2 ) - e N ( t s ) ) x ( t + Δt ) = x ( t s ) + α ( x ( t s + Δt 2 ) - x ( t s ) ) i N ( t + Δt ) = i N ( t s ) + α ( i N ( t s + Δt 2 ) - i N ( t s ) ) , 使与原来的仿真时标保持同步；
方法二：对混合框架进行由
Figure BDA00003714836400000413
时刻到
Figure BDA00003714836400000414
时刻的梯形法积分求解，即由
Figure BDA00003714836400000415
时刻对状态变量节点分析混合框架的计算方程Gglobalu=iglobal应用梯形法积分至
Figure BDA00003714836400000416
时刻，并计算得到
Figure BDA00003714836400000417
Figure BDA00003714836400000418
和
Figure BDA00003714836400000419
其中Ws和iN_hist参与计算方程Gglobalu=iglobal的求解；继续进行一步线性插值得到t+Δt时刻的值，此时
Figure BDA0000371483640000051
即计算 e N ( t + Δt ) = e N ( t s + Δt 2 ) + α ( e N ( t s + 3 Δt 2 ) - e N ( t s + Δt 2 ) ) x ( t + Δt ) = x ( t s + Δt 2 ) + α ( x ( t s + 3 Δt 2 ) - x ( t s + Δt 2 ) ) i N ( t + Δt ) = i N ( t s + Δt 2 ) + α ( i N ( t s + 3 Δt 2 ) - i N ( t s + Δt 2 ) ) , 使与原来的仿真时标保持同步。
2.根据权利要求1所述的面向状态变量节点分析混合框架的电磁暂态仿真插值方法，其特征在于，步骤2）所述的混合框架进行由t时刻到t+Δt时刻的梯形法积分求解，是由t时刻对状态变量节点分析混合框架的计算方程Gglobalu=iglobal应用梯形法积分至t+Δt时刻，并计算得到eN(t+Δt)、x(t+Δt)和iN(t+Δt)，其中由状态方程运用梯形法差分化得到的等效节点电导矩阵Ws和等效节点注入电流源列向量iN_hist参与计算方程Gglobalu=iglobal的求解。
3.根据权利要求1所述的面向状态变量节点分析混合框架的电磁暂态仿真插值方法，其特征在于，步骤4）所述的对混合框架进行由ts时刻到
Figure BDA0000371483640000053
时刻的半步长后向欧拉法积分求解，是由ts时刻对状态变量节点分析混合框架的计算方程Gglobalu=iglobal应用半步长后向欧拉法积分至
Figure BDA0000371483640000054
时刻，并计算得到
Figure BDA0000371483640000056
和
Figure BDA0000371483640000057
其中由状态方程运用半步长后向欧拉法差分化得到的等效节点电导矩阵W′s和等效节点注入电流源列向量i′N_hist参与计算方程Gglobalu=iglobal的求解。
4.根据权利要求1所述的面向状态变量节点分析混合框架的电磁暂态仿真插值方法，其特征在于，步骤7）所述的混合框架由时刻到
Figure BDA0000371483640000059
时刻采用两步半步长后向欧拉法进行重新初始化，是在
Figure BDA00003714836400000510
时刻对节点方程和状态方程采用半步长后向欧拉法积分至ts，得到eN(ts)、x(ts)和iN(ts)，其中W′s和i′N_hist参与计算方程Gglobalu=iglobal的求解；在ts时刻在采用半步长后向欧拉法积分至
Figure BDA00003714836400000511
得到
Figure BDA00003714836400000512
Figure BDA00003714836400000513
和
Figure BDA00003714836400000514
其中Ws′和i′N_hist参与计算方程Gglobalu=iglobal的求解，此时得到节点方程与状态方程中各个变量在时刻正确的值。
本发明的面向状态变量节点分析混合框架的电磁暂态仿真插值方法，实现了开关动作时刻系统状态的计算与重初始化，消除了开关动作之后数值振荡问题，同时通过合理的计算时序设计考虑了可能出现在同一步长内的多重开关动作问题，从而使混合仿真框架具有了精确处理大量电力电子开关快速动作的能力，极大地提升了状态变量节点分析混合仿真框架的仿真计算能力。本发明在电力系统电磁暂态仿真及含分布式电源和储能装置的智能配电网、局部孤岛运行电网的暂态仿真的实用计算方面具有重要意义。
附图说明
图1是节点方程与状态方程接口示意图；
图2是电感与电阻串联电路示意图；
图3是电容与电阻并联电路示意图；
图4是状态方程经梯形法差分后的等效电路示意图；
图5是本发明的面向状态变量节点分析混合框架的电磁暂态仿真插值方法的流程图；
图6是考虑光伏发电系统接入的欧盟低压微网算例结构图；
图7是光伏发电系统并网变压器高压侧A相电压波形；
图8是光伏发电系统并网变压器高压侧A相电流波形。
具体实施方式
下面结合实施例和附图对本发明的面向状态变量节点分析混合框架的电磁暂态仿真插值方法做出详细说明。
本发明的面向状态变量节点分析混合框架的电磁暂态仿真插值方法，属于电力系统仿真领域。TSDG（Transient Simulator for Distributed Generation Systems and Microgrid）是一种面向分布式发电微网系统暂态仿真的计算程序，实现了状态变量节点分析混合框架的电磁暂态仿真插值方法，本发明以欧盟低压微网算例作为测试算例，欧盟低压微网系统结构如图6所示。欧盟低压微网算例是在欧盟第五框架计划支持下的微网研究项目“Microgrids”提出的一个用于微网设计、仿真与测试的低压微网算例，系统中含有多种线路与负荷类型，可接入多种形式的分布式电源，充分体现了微网结构与运行的复杂性。作为测试，只接入了单级光伏发电系统。算例仿真时间为0.5s，仿真步长为2.5us，0.40s时光伏发电系统并网变压器高压侧A相接地短路，下面以欧盟低压微网算例作为测试算例详细说明本发明的一种面向状态变量节点分析混合框架的电磁暂态仿真插值方法，如图5所示，包括如下步骤：
1）基于节点方程和状态方程分别建立研究系统的电磁暂态仿真模型，本实施例是基于状态方程建立欧盟低压微网算例网络部分的电磁暂态仿真模型，基于节点方程建立欧盟低压微网算例光伏发电系统的电磁暂态仿真模型，确定节点方程与状态方程的接口变量，以联系节点电压eN为状态方程的输入量和以状态方程注入联系节点电流iN为状态方程的输出量，即u=eN，y=iN；
2）对由节点方程和状态方程构成的混合框架进行由t时刻到t+Δt时刻的梯形法积分求解；所述的混合框架进行由t时刻到t+Δt时刻的梯形法积分求解，是由t时刻对状态变量节点分析混合框架的计算方程Gglobalu=iglobal应用梯形法积分至t+Δt时刻，并计算得到eN(t+Δt)、x(t+Δt)和iN(t+Δt)，其中由状态方程运用梯形法差分化得到的等效节点电导矩阵Ws和等效节点注入电流源列向量iN_hist参与计算方程Gglobalu=iglobal的求解。
3）检测所有开关动作情况，如果不存在开关动作则转步骤2），否则通过线性插值的方法，由t和t+Δt时刻节点方程与状态方程中各值内插到最先动作开关ts时刻值，即计算
e N ( t s ) = e N ( t ) + α ( e N ( t + Δt ) - e N ( t ) ) x ( t s ) = x ( t ) + α ( x ( t + Δt ) - x ( t ) ) i N ( t s ) = i N ( t ) + α ( i N ( t + Δt ) - i N ( t ) ) ,
其中插值系数由开关电流最先反向的开关在t和t+Δt时刻电流值is(t)和is(t+Δt)确定，即 α = i s ( t ) i s ( t ) - i s ( t + Δt ) ;
4）改变开关状态，重新形成暂态计算矩阵，对混合框架进行由ts时刻到
Figure BDA0000371483640000073
时刻的半步长后向欧拉法积分求解；所述的对混合框架进行由ts时刻到
Figure BDA0000371483640000074
时刻的半步长后向欧拉法积分求解，是由ts时刻对状态变量节点分析混合框架的计算方程Gglobalu=iglobal应用半步长后向欧拉法积分至
Figure BDA0000371483640000075
时刻，并计算得到
Figure BDA0000371483640000077
和
Figure BDA0000371483640000078
其中由状态方程运用半步长后向欧拉法差分化得到的等效节点电导矩阵Ws′和等效节点注入电流源列向量i′N_hist参与计算方程Gglobalu=iglobal的求解。
5）重新检测所有开关动作情况，如果没有其他开关动作，则转入下一步，否则返回步骤4）；
6）运用线性插值方法，对混合框架中的节点方程和状态方程ts和
Figure BDA0000371483640000079
时刻值外插到
Figure BDA00003714836400000710
时刻值，此时α=-1.0，即计算 e N ( t s - Δt 2 ) = e N ( t s ) + α ( e N ( t s + Δt 2 ) - e N ( t s ) ) x ( t s - Δt 2 ) = x ( t s ) + α ( x ( t s + Δt 2 ) - x ( t s ) ) i N ( t s - Δt 2 ) = i N ( t s ) + α ( i N ( t s + Δt 2 ) - i N ( t s ) ) ;
7）对混合框架进行由
Figure BDA00003714836400000712
时刻到
Figure BDA00003714836400000713
时刻采用两步半步长后向欧拉法进行重新初始化；所述的混合框架由时刻到时刻采用两步半步长后向欧拉法进行重新初始化，是在
Figure BDA00003714836400000716
时刻对节点方程和状态方程采用半步长后向欧拉法积分至ts，得到eN(ts)、x(ts)和iN(ts)，其中Ws′和i′N_hist参与计算方程Gglobalu=iglobal的求解；在ts时刻在采用半步长后向欧拉法积分至
Figure BDA00003714836400000717
得到
Figure BDA00003714836400000718
Figure BDA00003714836400000719
和
Figure BDA00003714836400000720
其中Ws′和i′N_hist参与计算方程Gglobalu=iglobal的求解，此时得到节点方程与状态方程中各个变量在
Figure BDA00003714836400000721
时刻正确的值。
8）比较时刻
Figure BDA00003714836400000722
与时刻t+Δt的前后关系，若则采用如下方法一，否则采用如下方法二；
方法一：进行一步线性插值得到t+Δt时刻的值，此时即计算 e N ( t + Δt ) = e N ( t s ) + α ( e N ( t s + Δt 2 ) - e N ( t s ) ) x ( t + Δt ) = x ( t s ) + α ( x ( t s + Δt 2 ) - x ( t s ) ) i N ( t + Δt ) = i N ( i s ) + α ( i N ( t s + Δt 2 ) - i N ( t s ) ) , 使与原来的仿真时标保持同步；
方法二：对混合框架进行由
Figure BDA0000371483640000083
时刻到
Figure BDA0000371483640000084
时刻的梯形法积分求解，即由时刻对状态变量节点分析混合框架的计算方程Gglobalu=iglobal应用梯形法积分至
Figure BDA0000371483640000086
时刻，并计算得到
Figure BDA0000371483640000087
Figure BDA0000371483640000088
和
Figure BDA0000371483640000089
其中Ws和iN_hist参与计算方程Gglobalu=iglobal的求解；继续进行一步线性插值得到t+△t时刻的值，此时
Figure BDA00003714836400000810
即计算 e N ( t + Δt ) = e N ( t s + Δt 2 ) + α ( e N ( t s + 3 Δt 2 ) - e N ( t s + Δt 2 ) ) x ( t + Δt ) = x ( t s + Δt 2 ) + α ( x ( t s + 3 Δt 2 ) - x ( t s + Δt 2 ) ) i N ( t + Δt ) = i N ( t s + Δt 2 ) + α ( i N ( t s + 3 Δt 2 ) - i N ( t s + Δt 2 ) ) , 使与原来的仿真时标保持同步。
执行仿真计算的计算机硬件环境为Intel Core2Q84002.66GHz CPU，内存容量2GB；软件环境为Windows7操作系统。
附图7和图8比较了采用本发明面向状态变量节点分析混合框架的电磁暂态仿真插值方法的TSDG与商业仿真软件MATLAB/SimPowerSystems的仿真结果，为了便于观察与绘图，图中仅给出了0.30s-0.50s的仿真结果。从图中可以看出，MATLAB/SimPowerSystems仿真结果与TSDG的仿真结果在稳态与暂态过程中都能够完全吻合，二者的动态响应特性保持了高度一致，体现出了良好的仿真精度，充分验证了本发明提出的面向状态变量节点分析混合框架的电磁暂态仿真插值方法的正确性与可行性。
以上算例测试结果证明，本发明提出的面向状态变量节点分析混合框架的电磁暂态仿真插值方法具有良好的可行性与适用性，实现了开关动作时刻系统状态的计算与重初始化，消除了开关动作之后数值振荡问题，同时通过合理的计算时序设计考虑了可能出现在同一步长内的多重开关动作问题，从而使混合仿真框架具有了精确处理大量电力电子开关快速动作的能力，极大地提升了状态变量节点分析混合仿真框架的仿真计算能力。

Claims (4)

1.一种面向状态变量节点分析混合框架的电磁暂态仿真插值方法，其特征在于，包括如下步骤：
1）基于节点方程和状态方程分别建立研究系统的电磁暂态仿真模型，确定节点方程与状态方程的接口变量，以联系节点电压eN为状态方程的输入量和以状态方程注入联系节点电流iN为状态方程的输出量，即u=eN，y=iN；
2）对由节点方程和状态方程构成的混合框架进行由t时刻到t+Δt时刻的梯形法积分求解；
3）检测所有开关动作情况，如果不存在开关动作则转步骤2），否则通过线性插值的方法，由t和t+Δt时刻节点方程与状态方程中各值内插到最先动作开关ts时刻值，即计算
e N ( t s ) = e N ( t ) + α ( e N ( t + Δt ) - e N ( t ) ) x ( t s ) = x ( t ) + α ( x ( t + Δt ) - x ( t ) ) i N ( t s ) = i N ( t ) + α ( i N ( t + Δt ) - i N ( t ) ) ,
其中插值系数由开关电流最先反向的开关在t和t+Δt时刻电流值is(t)和is(t+Δt)确定，即 α = i s ( t ) i s ( t ) - i s ( t + Δt ) ;
4）改变开关状态，重新形成暂态计算矩阵，对混合框架进行由ts时刻到时刻的半步长后向欧拉法积分求解；
5）重新检测所有开关动作情况，如果没有其他开关动作，则转入下一步，否则返回步骤4）；
6）运用线性插值方法，对混合框架中的节点方程和状态方程ts和时刻值外插到
Figure FDA0000371483630000015
时刻值，此时α=-1.0，即计算 e N ( t s - Δt 2 ) = e N ( t s ) + α ( e N ( t s + Δt 2 ) - e N ( t s ) ) x ( t s - Δt 2 ) = x ( t s ) + α ( x ( t s + Δt 2 ) - x ( t s ) ) i N ( t s - Δt 2 ) = i N ( t s ) + α ( i N ( t s + Δt 2 ) - i N ( t s ) ) ;
7）对混合框架进行由
Figure FDA0000371483630000017
时刻到
Figure FDA0000371483630000018
时刻采用两步半步长后向欧拉法进行重新初始化；
8）比较时刻
Figure FDA0000371483630000019
与时刻t+Δt的前后关系，若
Figure FDA00003714836300000110
则采用如下方法一，否则采用如下方法二；
方法一：进行一步线性插值得到t+Δt时刻的值，此时即计算 e N ( t + Δt ) = e N ( t s ) + α ( e N ( t s + Δt 2 ) - e N ( t s ) ) x ( t + Δt ) = x ( t s ) + α ( x ( t s + Δt 2 ) - x ( t s ) ) i N ( t + Δt ) = i N ( t s ) + α ( i N ( t s + Δt 2 ) - i N ( t s ) ) , 使与原来的仿真时标保持同步；
方法二：对混合框架进行由
Figure FDA0000371483630000023
时刻到
Figure FDA0000371483630000024
时刻的梯形法积分求解，即由
Figure FDA0000371483630000025
时刻对状态变量节点分析混合框架的计算方程Gglobalu=iglobal应用梯形法积分至
Figure FDA0000371483630000026
时刻，并计算得到
Figure FDA0000371483630000027
Figure FDA0000371483630000028
和
Figure FDA0000371483630000029
其中Ws和iN_hist参与计算方程Gglobalu=iglobal的求解；继续进行一步线性插值得到t+Δt时刻的值，此时
Figure FDA00003714836300000210
即计算 e N ( t + Δt ) = e N ( t s + Δt 2 ) + α ( e N ( t s + 3 Δt 2 ) - e N ( t s + Δt 2 ) ) x ( t + Δt ) = x ( t s + Δt 2 ) + α ( x ( t s + 3 Δt 2 ) - x ( t s + Δt 2 ) ) i N ( t + Δt ) = i N ( t s + Δt 2 ) + α ( i N ( t s + 3 Δt 2 ) - i N ( t s + Δt 2 ) ) , 使与原来的仿真时标保持同步。
2.根据权利要求1所述的面向状态变量节点分析混合框架的电磁暂态仿真插值方法，其特征在于，步骤2）所述的混合框架进行由t时刻到t+Δt时刻的梯形法积分求解，是由t时刻对状态变量节点分析混合框架的计算方程Gglobalu=iglobal应用梯形法积分至t+Δt时刻，并计算得到eN(t+Δt)、x(t+Δt)和iN(t+Δt)，其中由状态方程运用梯形法差分化得到的等效节点电导矩阵Ws和等效节点注入电流源列向量iN_hist参与计算方程Gglobalu=iglobal的求解。
3.根据权利要求1所述的面向状态变量节点分析混合框架的电磁暂态仿真插值方法，其特征在于，步骤4）所述的对混合框架进行由ts时刻到时刻的半步长后向欧拉法积分求解，是由ts时刻对状态变量节点分析混合框架的计算方程Gglobalu=iglobal应用半步长后向欧拉法积分至
Figure FDA00003714836300000213
时刻，并计算得到
Figure FDA00003714836300000214
Figure FDA00003714836300000215
和
Figure FDA00003714836300000216
其中由状态方程运用半步长后向欧拉法差分化得到的等效节点电导矩阵W′s和等效节点注入电流源列向量i′N_hist参与计算方程Gglobalu=iglobal的求解。
4.根据权利要求1所述的面向状态变量节点分析混合框架的电磁暂态仿真插值方法，其特征在于，步骤7）所述的混合框架由时刻到
Figure FDA0000371483630000032
时刻采用两步半步长后向欧拉法进行重新初始化，是在
Figure FDA0000371483630000033
时刻对节点方程和状态方程采用半步长后向欧拉法积分至ts，得到eN(ts)、x(ts)和iN(ts)，其中W′s和i′N_hist参与计算方程Gglobalu=iglobal的求解；在ts时刻在采用半步长后向欧拉法积分至
Figure FDA0000371483630000034
得到
Figure FDA0000371483630000035
Figure FDA0000371483630000036
和
Figure FDA0000371483630000037
其中W′s和i′N_hist参与计算方程Gglobalu=iglobal的求解，此时得到节点方程与状态方程中各个变量在
Figure FDA0000371483630000038
时刻正确的值。
CN201310374258.3A 2013-08-23 2013-08-23 面向状态变量节点分析混合框架的电磁暂态仿真插值方法 Active CN103455668B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201310374258.3A CN103455668B (zh) 2013-08-23 2013-08-23 面向状态变量节点分析混合框架的电磁暂态仿真插值方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201310374258.3A CN103455668B (zh) 2013-08-23 2013-08-23 面向状态变量节点分析混合框架的电磁暂态仿真插值方法

Publications (2)

Publication Number Publication Date
CN103455668A true CN103455668A (zh) 2013-12-18
CN103455668B CN103455668B (zh) 2016-01-13

Family

ID=49738024

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201310374258.3A Active CN103455668B (zh) 2013-08-23 2013-08-23 面向状态变量节点分析混合框架的电磁暂态仿真插值方法

Country Status (1)

Country Link
CN (1) CN103455668B (zh)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105989211A (zh) * 2015-03-01 2016-10-05 范圣韬 基于节点分析法的电压定义支路实现方法
CN106446428A (zh) * 2016-09-29 2017-02-22 全球能源互联网研究院 一种开关电路电磁暂态分析方法及分析装置
CN108563588A (zh) * 2018-03-18 2018-09-21 天津大学 基于fpga的有源配电网实时仿真器多速率接口设计方法
CN109508492A (zh) * 2018-11-10 2019-03-22 东北电力大学 一种交直流混杂模式下的y/δ变压器铁芯振动计算方法
CN112487629A (zh) * 2020-11-25 2021-03-12 南方电网科学研究院有限责任公司 考虑多重事件发生的电磁暂态仿真方法、装置以及设备

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102708250A (zh) * 2012-05-10 2012-10-03 天津大学 基于树形分层双向迭代的电力系统数字混合仿真方法
CN102819641A (zh) * 2012-08-08 2012-12-12 天津大学 适于电磁暂态仿真的大规模配电网络整体模型化简方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102708250A (zh) * 2012-05-10 2012-10-03 天津大学 基于树形分层双向迭代的电力系统数字混合仿真方法
CN102819641A (zh) * 2012-08-08 2012-12-12 天津大学 适于电磁暂态仿真的大规模配电网络整体模型化简方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
CHRISTIAN DUFOUR 等: "A Combined State-Space Nodal Method for the Simulation of Power System Transients", 《IEEE TRANSACTIONS ON POWER DELIVERY》, vol. 26, no. 2, 30 April 2011 (2011-04-30) *
赵帅 等: "一种考虑多重开关动作的变步长电磁暂态仿真算法研究", 《电工技术学报》, vol. 27, no. 3, 31 March 2012 (2012-03-31) *

Cited By (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105989211A (zh) * 2015-03-01 2016-10-05 范圣韬 基于节点分析法的电压定义支路实现方法
CN105989211B (zh) * 2015-03-01 2021-12-21 范圣韬 一种基于节点法实现电压定义电流控制支路元件的方法
CN106446428A (zh) * 2016-09-29 2017-02-22 全球能源互联网研究院 一种开关电路电磁暂态分析方法及分析装置
CN108563588A (zh) * 2018-03-18 2018-09-21 天津大学 基于fpga的有源配电网实时仿真器多速率接口设计方法
CN108563588B (zh) * 2018-03-18 2021-02-02 天津大学 基于fpga的有源配电网实时仿真器多速率接口设计方法
CN109508492A (zh) * 2018-11-10 2019-03-22 东北电力大学 一种交直流混杂模式下的y/δ变压器铁芯振动计算方法
CN112487629A (zh) * 2020-11-25 2021-03-12 南方电网科学研究院有限责任公司 考虑多重事件发生的电磁暂态仿真方法、装置以及设备
CN112487629B (zh) * 2020-11-25 2024-05-31 南方电网科学研究院有限责任公司 考虑多重事件发生的电磁暂态仿真方法、装置以及设备

Also Published As

Publication number Publication date
CN103455668B (zh) 2016-01-13

Similar Documents

Publication Publication Date Title
Liu et al. Online voltage stability assessment for load areas based on the holomorphic embedding method
Gurrala et al. Parareal in time for fast power system dynamic simulations
Dagbagi et al. ADC-based embedded real-time simulator of a power converter implemented in a low-cost FPGA: Application to a fault-tolerant control of a grid-connected voltage-source rectifier
CN103455668B (zh) 面向状态变量节点分析混合框架的电磁暂态仿真插值方法
CN103646152A (zh) 一种基于矩阵指数的电力系统电磁暂态仿真方法
Rakpenthai et al. State estimation of power system considering network parameter uncertainty based on parametric interval linear systems
CN104298809B (zh) 一种基于矩阵指数电磁暂态仿真的非线性建模求解方法
CN102184297B (zh) 适于微网暂态并行仿真的电气与控制系统解耦预测方法
CN103077268B (zh) 面向电力系统电磁暂态仿真的状态空间自动建模方法
CN104217074B (zh) 一种基于矩阵指数的电磁暂态隐式降阶仿真方法
CN103440374A (zh) 基于状态变量节点分析混合框架的电磁暂态仿真建模方法
CN103853052A (zh) 一种核电站反应堆控制系统的设计方法
CN104967114B (zh) 电网负荷实时数字建模方法及系统
Belaouar et al. An asymptotically stable semi-lagrangian scheme in the quasi-neutral limit
CN106844900B (zh) 电磁暂态仿真系统的搭设方法
CN104779613B (zh) 基于试验的含变流器电力元件等效建模方法
Chitturi et al. Comparing performance of Prony analysis and matrix pencil method for monitoring power system oscillations
CN105302979A (zh) 两相流体网络模型中阀门组的建模方法和系统
Abourida et al. Hardware-in-the-loop simulation of electric systems and power electronics on FPGA using physical modeling
Hayes et al. Viable computation of the largest Lyapunov characteristic exponent for power systems
Zhang et al. An implicit algorithm for high-order DG/FV schemes for compressible flows on 2D arbitrary grids
CN108090846B (zh) 一种电网低频振荡案例库的构建方法及装置
Abhyankar et al. Parallel-in-space-and-time scheme for implicitly coupled electromechanical and electromagnetic transients simulation
Wang et al. Uncertainty importance measure of individual components in multi-state systems
CN103779863B (zh) 一种在极坐标系下减小区间潮流保守性的方法

Legal Events

Date Code Title Description
C06 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
TR01 Transfer of patent right

Effective date of registration: 20210827

Address after: 510700 3rd, 4th and 5th floors of building J1 and 3rd floor of building J3, No.11 Kexiang Road, Science City, Luogang District, Guangzhou City, Guangdong Province

Patentee after: China Southern Power Grid Research Institute Co.,Ltd.

Patentee after: Tianjin University

Address before: West Tower, Yuedian building, No. 6 and 8, shuijungang, Dongfeng East Road, Yuexiu District, Guangzhou, Guangdong

Patentee before: CSG POWER GRID TECHNOLOGY RESEARCH CENTER

Patentee before: Tianjin University

TR01 Transfer of patent right