CN113221298A - 一种机电暂态过程的仿真方法及系统 - Google Patents

一种机电暂态过程的仿真方法及系统 Download PDF

Info

Publication number
CN113221298A
CN113221298A CN202110434097.7A CN202110434097A CN113221298A CN 113221298 A CN113221298 A CN 113221298A CN 202110434097 A CN202110434097 A CN 202110434097A CN 113221298 A CN113221298 A CN 113221298A
Authority
CN
China
Prior art keywords
differential
equation
simulation
algebraic
power system
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
CN202110434097.7A
Other languages
English (en)
Other versions
CN113221298B (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.)
China Southern Power Grid Co Ltd
Research Institute of Southern Power Grid Co Ltd
Original Assignee
China Southern Power Grid Co Ltd
Research Institute of 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 China Southern Power Grid Co Ltd, Research Institute of Southern Power Grid Co Ltd filed Critical China Southern Power Grid Co Ltd
Priority to CN202110434097.7A priority Critical patent/CN113221298B/zh
Publication of CN113221298A publication Critical patent/CN113221298A/zh
Application granted granted Critical
Publication of CN113221298B publication Critical patent/CN113221298B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/10—Geometric CAD
    • G06F30/18—Network design, e.g. design based on topological or interconnect aspects of utility systems, piping, heating ventilation air conditioning [HVAC] or cabling
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/13—Differential equations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02E—REDUCTION OF GREENHOUSE GAS [GHG] EMISSIONS, RELATED TO ENERGY GENERATION, TRANSMISSION OR DISTRIBUTION
    • Y02E60/00—Enabling technologies; Technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Geometry (AREA)
  • Evolutionary Computation (AREA)
  • Computer Hardware Design (AREA)
  • Data Mining & Analysis (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Operations Research (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Supply And Distribution Of Alternating Current (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明公开了一种机电暂态过程的仿真方法及系统，方法利用所述电力系统初始化时的运行参数，获得初始仿真结果；基于牛顿积分算法，根据所述初始仿真结果计算所述第一微分‑代数方程组的雅可比矩阵，并在相邻两次迭代的雅可比矩阵解的残差小于第一预设阈值时，不更新雅可比矩阵；基于牛顿积分算法，对所述第二微分‑代数方程组进行迭代求解，在相邻两次迭代的状态变量的残差小于第二预设阈值，且相邻两次迭代的网络变量的残差小于第三预设阈值时，计算收敛，进行下一步仿真计算。本发明通过采用牛顿积分算法求解微分‑代数方程组，确保了仿真的精度，在迭代过程中不更新雅克比矩阵，能够缩短机电暂态过程的仿真计算时间。

Description

一种机电暂态过程的仿真方法及系统
技术领域
本发明涉及电子系统仿真领域，尤其涉及一种机电暂态过程的仿真方法及系统。
背景技术
现代电力系统已进入大系统、超高压、跨区域联网和远距离输电的新时期，对电力系统进行详细与准确的建模、实施高效和快速的仿真是一项确保电力系统安全运行的重要手段。
大部分电力系统暂态仿真问题都是刚性的，所以机电暂态仿真程序多采用隐式积分法，它的关键是牛顿积分算法。经研究，在暂态仿真过程中，牛顿积分法每一步长都要更新雅可比矩阵，由于雅可比矩阵的频繁更新必会引起矩阵的重新分解，而且这个矩阵并不是稀疏矩阵，导致仿真计算时间长。
发明内容
本发明提供一种机电暂态过程的仿真方法及系统，以解决现有技术仿真计算时间长的问题，本发明能够缩短机电暂态过程的仿真计算时间，满足实时监控电力系统的需求。
为解决以上的技术问题，本发明实施例提供了一种机电暂态过程的仿真方法，包括：
利用电力系统的元件参数和网络拓扑结构，形成描述所述电力系统机电暂态过程的第一微分-代数方程组；
利用所述电力系统初始化时的运行参数，获得初始仿真结果，其中，所述初始仿真结果包括状态变量的初始值和网络变量的初始值；
基于牛顿积分算法，根据所述初始仿真结果计算所述第一微分-代数方程组的雅可比矩阵，并在相邻两次迭代的雅可比矩阵解的残差小于第一预设阈值时，不更新雅可比矩阵，获得所述雅可比矩阵、所述第一微分-代数方程组的第一状态变量和所述第一微分-代数方程组的第一网络变量；
利用所有设备的注入电流，形成描述所述电力系统机电暂态过程的第二微分-代数方程组；
基于牛顿积分算法，根据所述雅可比矩阵、所述第一状态变量和所述第一网络变量对所述第二微分-代数方程组进行迭代求解，在相邻两次迭代的状态变量的残差小于第二预设阈值，且相邻两次迭代的网络变量的残差小于第三预设阈值时，获得所述第二微分-代数方程组的第二状态变量和所述第二微分-代数方程组的第二网络变量；
当仿真时间达到仿真终止时间，输出仿真结果；其中，所述仿真结果包括所述电力系统各个节点的电压、功率、发电机功角曲线和发电机相对摇摆角曲线。
进一步地，所述第一微分-代数方程组如公式(1)所示：
Figure BDA0003031259990000021
其中，微分方程表示电力系统元件的动态特性，是电力系统的状态方程；代数方程表示电力系统元件的静态特性，主要是系统的网络方程；x为n个状态向量，V为m个代数向量，微分代数方程组的阶数为n和m之和；Y是电力系统的导纳矩阵。
进一步地，所述利用所述电力系统初始化时的运行参数，获得初始仿真结果，其中，所述初始仿真结果包括状态变量的初始值和网络变量的初始值，具体为：
利用所述电力系统初始化时的运行参数，采用简单欧拉法对初始状态变量进行预测，获得状态变量的初始值，则所述状态变量的初始值
Figure BDA0003031259990000022
为：
Figure BDA0003031259990000023
其中，h为仿真步长，xn为第n步的状态变量，Vn为第n步的网络变量；
利用所述电力系统初始化时的运行参数，采用几何预测法对初始网络变量进行预测，获得网络变量的初始值，则所述网络变量的初始值
Figure BDA0003031259990000024
为：
Figure BDA0003031259990000025
其中，Vn-1为第n-1步的网络变量，Vn为第n步的网络变量。
进一步地，，所述基于牛顿积分算法，根据所述初始仿真结果计算所述第一微分-代数方程组的雅可比矩阵，具体为：
对所述第一微分-代数方程组采用隐式梯形积分法进行差分，得到公式(2)和(3)：
Figure BDA0003031259990000031
YVn+1＝I(xn+1，Vn+1) (3)
其中，Y是电力系统的导纳矩阵，h是步长；
令公式(2)改成为公式(4)：
Figure BDA0003031259990000032
根据所述状态变量的初始值和所述网络变量的初始值，采用牛顿积分算法求解F(xn+1，Vn+1)＝0，计算所述第一微分-代数方程组的雅可比矩阵。
进一步地，所述利用所有设备的注入电流，形成描述所述电力系统机电暂态过程的第二微分-代数方程组，具体为：
对公式(2)和公式(3)采用牛顿积分算法求解，得到公式(10)和(14)：
Figure BDA0003031259990000033
Figure BDA0003031259990000034
其中，AG、BG、CG和YG是雅克比矩阵，雅可比矩阵AG、BG、CG和YG是(x，V)的函数；
利用所有设备的注入电流，得到公式(18)：
Figure BDA0003031259990000035
将公式(18)代替公式(14)，形成所述第二微分-代数方程组，所述第二微分-代数方程组如公式(10)和公式(18)所示：
Figure BDA0003031259990000036
Figure BDA0003031259990000037
其中，YGn是电力系统所有设备对应于非凸极效应的常数矩阵。
进一步地，所述相邻两次迭代的状态变量的残差具体为：
Figure BDA0003031259990000041
其中，
Figure BDA0003031259990000042
为第k+1次迭代的状态变量，
Figure BDA0003031259990000043
为第k次迭代的状态变量；
所述相邻两次迭代的网络变量的残差具体为：
Figure BDA0003031259990000044
其中，
Figure BDA0003031259990000045
为第k+1次迭代的状态变量，
Figure BDA0003031259990000046
为第k次迭代的状态变量。
所述方法还包括：
当仿真时间未达到仿真终止时间，继续进行下一步长仿真计算。
相应地，本发明实施例还提供了一种机电暂态过程的仿真系统，包括：
第一仿真模型单元括状态变量的初始值和网件参数和网络拓扑结构，形成描述所述电力系统机电暂态过程的第一微分-代数方程组；
初始化单元：用于利用所述电力系统初始化时的运行参数，获得初始仿真结果，其中，所述初始仿真结果包括状态变量的初始值和网络变量的初始值；
雅克比矩阵单元：用于基于牛顿积分算法，根据所述初始仿真结果计算所述第一微分-代数方程组的雅可比矩阵，并在相邻两次迭代的雅可比矩阵解的残差小于第一预设阈值时，不更新雅可比矩阵，获得所述雅可比矩阵、所述第一微分-代数方程组的第一状态变量和所述第一微分-代数方程组的第一网络变量；
第二仿真模型单元：用于利用所有设备的注入电流，形成描述所述电力系统机电暂态过程的第二微分-代数方程组；
求解单元：用于基于牛顿积分算法，根据所述雅可比矩阵、所述第一状态变量和所述第一网络变量对所述第二微分-代数方程组进行迭代求解，在相邻两次迭代的状态变量的残差小于第二预设阈值，且相邻两次迭代的网络变量的残差小于第三预设阈值时，获得所述第二微分-代数方程组的第二状态变量和所述第二微分-代数方程组的第二网络变量；
输出单元：用于当仿真时间达到仿真终止时间，输出仿真结果；其中，所述仿真结果包括所述电力系统各个节点的电压、功率、发电机功角曲线和发电机相对摇摆角曲线。
进一步地，所述第一微分-代数方程组如公式(1)所示：
Figure BDA0003031259990000051
其中，微分方程表示电力系统元件的动态特性，是电力系统的状态方程；代数方程表示电力系统元件的静态特性，主要是系统的网络方程；x为n个状态向量，V为m个代数向量，微分代数方程组的阶数为n和m之和；Y是电力系统的导纳矩阵。
进一步地，所述相邻两次迭代的状态变量的残差具体为：
Figure BDA0003031259990000052
其中，
Figure BDA0003031259990000053
为第k+1次迭代的状态变量，
Figure BDA0003031259990000054
为第k次迭代的状态变量；
所述相邻两次迭代的网络变量的残差具体为：
Figure BDA0003031259990000055
其中，
Figure BDA0003031259990000056
为第k+1次迭代的状态变量，
Figure BDA0003031259990000057
为第k次迭代的状态变量。
本发明实施例提供的一种机电暂态过程的仿真方法及系统，通过对状态变量的初始值和网络变量的初始值进行预测，能精确获取状态变量的初始值和网络变量的初始值，可以有效降低一个时间步长中的迭代次数，缩短仿真计算时间；在相邻两次迭代的雅可比矩阵解的残差小于第一预设阈值时，不更新雅可比矩阵，减少了更新雅克比矩阵的占用时间，缩短仿真计算时间；本发明实施例通过采用牛顿积分算法求解微分-代数方程组，确保了仿真的精度，在迭代过程中不更新雅克比矩阵，能够缩短机电暂态过程的仿真计算时间，满足实时监控电力系统的需求。
附图说明
图1是本发明实施例提供的一种机电暂态过程的仿真方法流程示意图；
图2是本发明实施例提供的最大发电机相对转子角图；
图3是本发明实施例提供的一种机电暂态过程的仿真系统。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
电力系统的各种元件都可以由数学模型表示出来，其暂态过程可用如下形式的微分代数方程组描述：
Figure BDA0003031259990000061
其中，微分方程表示电力系统元件的动态特性，是系统的状态方程；代数方程表示电力系统元件的静态特性，主要是系统的网络方程；x为n个状态向量(微分变量)；V为m个代数向量(代数变量)。微分代数方程组的阶数为n和m之和。
对于公式(1)所描述的电力系统微分代数方程组，一般可以采用隐式梯形积分方法进行微分方程的求解，假设初始预测值为
Figure BDA0003031259990000062
采用隐式梯形积分法进行差分后，公式(1)变为求解以下两个方程：
Figure BDA0003031259990000063
YVn+1＝I(xn+1，Vn+1) (3)
其中Y是系统的导纳矩阵，h是步长；
公式(2)和(3)可以改写成：
Figure BDA0003031259990000064
G(xn+1，Vn+1)＝YVn+1-I(xn+1，Vn+1) (5)
用牛顿法的迭代步骤来求解公式(4)和(5)，将公式(4)和(5)重写为：
F(xn+1，Vn+1)＝0 (6)
G(xn+1，Vn+1)＝0 (7)
对公式(6)和(7)用牛顿法重写为：
Figure BDA0003031259990000065
Figure BDA0003031259990000066
从公式(8)得：
Figure BDA0003031259990000067
将等式(10)代入公式(9)得：
Figure BDA0003031259990000071
或；
Figure BDA0003031259990000072
从公式(5)得：
Figure BDA0003031259990000073
因此，等式(12)可以降为：
Figure BDA0003031259990000074
注意等式(14)右边只包含已知的变量，因此，母线电压可从等式(14)直接求解。母线电压已知后，状态变量可从等式(10)求得。当残余向量达到误差范围或等式(14)右端(等值电流注入向量)在误差范围之内几乎不变，迭代结束。
从上述讨论可以看出，精确牛顿法比较慢，原因在于：雅可比矩阵AG，BG，CG和YG是(x，V)的函数，每一次迭代都必须更新，这样对大系统会很慢；网络矩阵(等式(14)的左侧)是雅可比矩阵的函数，因此每次迭代时，都必须重新分解。
针对以上问题，本发明实施例提供一种机电暂态过程的仿真方法，可以使雅可比矩阵在一个时间步长或几个时间步长中，保持不变，在保证精度的同时，能够缩短机电暂态过程的仿真计算时间。以下将结合图1，对本发明实施例作进一步说明。
参见图1，是本发明实施例提供的一种机电暂态过程的仿真方法流程示意图，所述方法包括：
S1，利用电力系统的元件参数和网络拓扑结构，形成描述所述电力系统机电暂态过程的第一微分-代数方程组；
具体的，步骤S1所述的第一微分-代数方程组如公式(1)所示：
Figure BDA0003031259990000075
其中，微分方程表示电力系统元件的动态特性，是电力系统的状态方程；代数方程表示电力系统元件的静态特性，主要是系统的网络方程；x为n个状态向量，V为m个代数向量，微分代数方程组的阶数为n和m之和；Y是系统的导纳矩阵。
S2，利用所述电力系统初始化时的运行参数，获得初始仿真结果，其中，所述初始仿真结果包括状态变量的初始值和网络变量的初始值；
具体的，所述利用所述电力系统初始化时的运行参数，获得初始仿真结果，其中，所述初始仿真结果包括状态变量的初始值和网络变量的初始值，具体为：
利用所述电力系统初始化时的运行参数，采用简单欧拉法对初始状态变量进行预测，获得状态变量的初始值，则所述状态变量的初始值
Figure BDA0003031259990000081
为：
Figure BDA0003031259990000082
其中，h为仿真步长，xn为第n步的状态变量，Vn为第n步的网络变量；
利用所述电力系统初始化时的运行参数，采用几何预测法对初始网络变量进行预测，获得网络变量的初始值，则所述网络变量的初始值
Figure BDA0003031259990000083
为：
Figure BDA0003031259990000084
其中，Vn-1为第n-1步的网络变量，Vn为第n步的网络变量。
在一个具体的实施例中，在任意仿真时刻，精确预测状态变量的初始值和网络变量的初始值，可以有效降低一个时间步长中的迭代次数。将初始预测值代入微分-代数方程组求解，当初始预测值越接近最终解，所需的迭代次数越少，仿真计算时间越短。
S3，基于牛顿积分算法，根据所述初始仿真结果计算所述第一微分-代数方程组的雅可比矩阵，并在相邻两次迭代的雅可比矩阵解的残差小于第一预设阈值时，不更新雅可比矩阵，获得所述雅可比矩阵、所述第一微分-代数方程组的第一状态变量和所述第一微分-代数方程组的第一网络变量；
具体的，所述基于牛顿积分算法，根据所述初始仿真结果计算所述第一微分-代数方程组的雅可比矩阵，具体为：
对所述第一微分-代数方程组采用隐式梯形积分法进行差分，得到公式(2)和(3)：
Figure BDA0003031259990000085
YVn+1＝I(xn+1，Vn+1) (3)
其中，Y是电力系统的导纳矩阵，h是步长；
令公式(2)改成为公式(4)：
Figure BDA0003031259990000091
根据所述状态变量的初始值和所述网络变量的初始值，采用牛顿积分算法求解F(xn+1，Vn+1)＝0，计算所述第一微分-代数方程组的雅可比矩阵。
优选的，在相邻两次迭代的雅可比矩阵解的残差不小于第一预设阈值时，更新雅可比矩阵，并对公式(10)进行求解
S4，利用所有设备的注入电流，形成描述所述电力系统机电暂态过程的第二微分-代数方程组；
具体的，所述利用所有设备的注入电流，形成描述所述电力系统机电暂态过程的第二微分-代数方程组，具体为：
对公式(2)和公式(3)采用牛顿积分算法求解，得到公式(10)和(14)：
Figure BDA0003031259990000092
Figure BDA0003031259990000093
其中，AG、BG、CG和YG是雅克比矩阵，雅可比矩阵AG、BG、CG和YG是(x，V)的函数；
利用所有设备的注入电流，得到公式(18)：
Figure BDA0003031259990000094
将公式(18)代替公式(14)，形成所述第二微分-代数方程组，所述第二微分-代数方程组如公式(10)和公式(18)所示：
Figure BDA0003031259990000095
Figure BDA0003031259990000096
其中，YGn是电力系统所有设备对应于非凸极效应的常数矩阵。
需要说明的是，公式(10)和公式(18)通过结合上述公式(1)～公式(14)得出，以下结合公式(1)～公式(14)对本发明实施例作进一步说明。
对于公式(3)，发电机电流可写为：
ig(xg，vg)＝-Ygvg+igi(xg) (15)
其中xg，vg分别是发电机的状态矢量和机端电压。Yg是发电机的导纳阵。只有发电机考虑凸极效应时，它是xg的函数。总的来说，可以将Yg分解成两部分：一个对应于非凸极效应的常数矩阵Ygn，另一个包含所有凸极效应的矩阵Ygs。
这样公式(15)变成：
ig(xg，vg)＝-Ygnvg-Ygsvg+igi(xg) (16)
对所有其它设备，可以写出类似于等式(15)。因此，所有的电流注入可写为
I(x，V)＝-YGV+Ii(x)＝-YGnV-YGsV+Ii(x) (17)
从公式(15)和(17)，
YG＝YGn+YGs
Figure BDA0003031259990000101
让
Figure BDA0003031259990000102
Figure BDA0003031259990000103
公式(13)可写成：
Figure BDA0003031259990000104
在一个具体实施例中，采用牛顿积分算法联合求解公式(10)和(18)，注意公式(18)引入了I2(k+1)变量，它又和未知电压向量相关联。因此在计算中用前一次迭代的电压值来替代当前的电压值以加快速度。公式(18)可写为：
Figure BDA0003031259990000105
公式(19)是网络等式(18)的实际电流注入向量。
如果在一个时间步长中，雅克比矩阵解的残差小于第一预设阈值，不更新雅可比矩阵，这是本发明实施例提供一种机电暂态过程的仿真方法的关键。因为更新雅可比矩阵会占用很多时间。
S5，基于牛顿积分算法，根据所述雅可比矩阵、所述第一状态变量和所述第一网络变量对所述第二微分-代数方程组进行迭代求解，在相邻两次迭代的状态变量的残差小于第二预设阈值，且相邻两次迭代的网络变量的残差小于第三预设阈值时，获得所述第二微分-代数方程组的第二状态变量和所述第二微分-代数方程组的第二网络变量；
具体的，所述相邻两次迭代的状态变量的残差具体为：
Figure BDA0003031259990000111
其中，
Figure BDA0003031259990000112
为第k+1次迭代的状态变量，
Figure BDA0003031259990000113
为第k次迭代的状态变量；
所述相邻两次迭代的网络变量的残差具体为：
Figure BDA0003031259990000114
其中，
Figure BDA0003031259990000115
为第k+1次迭代的状态变量，
Figure BDA0003031259990000116
为第k次迭代的状态变量。
S6，当仿真时间达到仿真终止时间，输出仿真结果；其中，所述仿真结果包括所述电力系统各个节点的电压、功率、发电机功角曲线和发电机相对摇摆角曲线。
优选的，所述方法还包括：当仿真时间未达到仿真终止时间，继续进行下一步长仿真计算。
本发明实施例提供的机电暂态过程的仿真方法可用于机电暂态仿真或中长期动态仿真。在一个具体实施例中，将本发明实施例提供的一种机电暂态过程的仿真方法应用于电力系统机电暂态仿真程序，将其与牛顿法的电力系统机电暂态仿真程序、西门子的商用电力系统机电暂态仿真程序即PSS/E分别进行仿真，其中，仿真系统由3862个节点，414个发电机组成。仿真时间为10秒，仿真步长0.02秒。其仿真的计算时间比较如下：
牛顿法 本发明实施例 PSS/E
CPU时间(秒) 25.6 8.1 13.5
加速比 3.1倍 1 1.6倍
在另一个具体实施例中，使用现有技术中的传统RK4算法、DIgSILENT、显式投影算法和隐式投影算法分别进行仿真，仿真结果如下：
Figure BDA0003031259990000117
Figure BDA0003031259990000121
可见，相对于现有技术的牛顿积分算法的电力系统机电暂态仿真，本发明实施例提供的机电暂态仿真能够缩短机电暂态过程的仿真计算时间。
参见图2，是本发明实施例提供的最大发电机相对转子角图；
在一个具体的实施例中，将本发明实施例提供的一种机电暂态过程的仿真方法应用于电力系统机电暂态仿真，将其与牛顿法的电力系统机电暂态仿真程序分别进行仿真，得到的最大发电机相对转子角如图2，可以看出两条曲线几乎重合，说明本发明实施例提供的一种机电暂态过程的仿真方法保证了计算精度。
由以上测试的结果，表明本发明实施例提供的一种机电暂态过程的仿真方法在保证仿真精度的同时，能够缩短机电暂态过程的仿真计算时间，为电力系统实时监控打下了坚实的基础。
参见图3，是本发明实施例提供的一种机电暂态过程的仿真系统，包括：
第一仿真模型单元11：用于利用电力系统的元件参数和网络拓扑结构，形成描述所述电力系统机电暂态过程的第一微分-代数方程组；
初始化单元12：用于利用所述电力系统初始化时的运行参数，获得初始仿真结果，其中，所述初始仿真结果包括状态变量的初始值和网络变量的初始值；
雅克比矩阵单元13：用于基于牛顿积分算法，根据所述初始仿真结果计算所述第一微分-代数方程组的雅可比矩阵，并在相邻两次迭代的雅可比矩阵解的残差小于第一预设阈值时，不更新雅可比矩阵，获得所述雅可比矩阵、所述第一微分-代数方程组的第一状态变量和所述第一微分-代数方程组的第一网络变量；
第二仿真模型单元14：用于利用所有设备的注入电流，形成描述所述电力系统机电暂态过程的第二微分-代数方程组；
求解单元15：用于基于牛顿积分算法，根据所述雅可比矩阵、所述第一状态变量和所述第一网络变量对所述第二微分-代数方程组进行迭代求解，在相邻两次迭代的状态变量的残差小于第二预设阈值，且相邻两次迭代的网络变量的残差小于第三预设阈值时，获得所述第二微分-代数方程组的第二状态变量和所述第二微分-代数方程组的第二网络变量；
输出单元16：用于当仿真时间达到仿真终止时间，输出仿真结果；其中，所述仿真结果包括所述电力系统各个节点的电压、功率、发电机功角曲线和发电机相对摇摆角曲线。
具体的，所述第一微分-代数方程组如公式(1)所示：
Figure BDA0003031259990000131
其中，微分方程表示电力系统元件的动态特性，是电力系统的状态方程；代数方程表示电力系统元件的静态特性，主要是系统的网络方程；x为n个状态向量，V为m个代数向量，微分代数方程组的阶数为n和m之和；Y是电力系统的导纳矩阵。
具体的，初始化单元12具体为：
用于利用所述电力系统初始化时的运行参数，采用简单欧拉法对初始状态变量进行预测，获得状态变量的初始值，则所述状态变量的初始值
Figure BDA0003031259990000132
为：
Figure BDA0003031259990000133
其中，h为仿真步长，xn为第n步的状态变量，Vn为第n步的网络变量；
用于利用所述电力系统初始化时的运行参数，采用几何预测法对初始网络变量进行预测，获得网络变量的初始值，则所述网络变量的初始值
Figure BDA0003031259990000141
为：
Figure BDA0003031259990000142
其中，Vn-1为第n-1步的网络变量，Vn为第n步的网络变量。
具体的，所述基于牛顿积分算法，根据所述初始仿真结果计算所述第一微分-代数方程组的雅可比矩阵，具体为：
对所述第一微分-代数方程组采用隐式梯形积分法进行差分，得到公式(2)和(3)：
Figure BDA0003031259990000143
YVn+1＝I(xn+1，Vn+1) (3)
其中，Y是电力系统的导纳矩阵，h是步长；
令公式(2)改成为公式(4)：
Figure BDA0003031259990000144
根据所述状态变量的初始值和所述网络变量的初始值，采用牛顿积分算法求解F(xn+1，Vn+1)＝0，计算所述第一微分-代数方程组的雅可比矩阵。
具体的，第二仿真模型单元14具体为：
用于对公式(2)和公式(3)采用牛顿积分算法求解，得到公式(10)和(14)：
Figure BDA0003031259990000145
Figure BDA0003031259990000146
其中，AG、BG、CG和YG是雅克比矩阵，雅可比矩阵AG、BG、CG和YG是(x，V)的函数；
用于利用所有设备的注入电流，得到公式(18)：
Figure BDA0003031259990000147
用于将公式(18)代替公式(14)，形成所述第二微分-代数方程组，所述第二微分-代数方程组如公式(10)和公式(18)所示：
Figure BDA0003031259990000148
Figure BDA0003031259990000151
其中，YGn是电力系统所有设备对应于非凸极效应的常数矩阵。
具体的，所述相邻两次迭代的状态变量的残差具体为：
Figure BDA0003031259990000152
其中，
Figure BDA0003031259990000153
为第k+1次迭代的状态变量，
Figure BDA0003031259990000154
为第k次迭代的状态变量；
所述相邻两次迭代的网络变量的残差具体为：
Figure BDA0003031259990000155
其中，
Figure BDA0003031259990000156
为第k+1次迭代的状态变量，
Figure BDA0003031259990000157
为第k次迭代的状态变量。
具体的，所述系统还包括：
用于当仿真时间未达到仿真终止时间，继续进行下一步长仿真计算。
本发明实施例提供的一种机电暂态过程的仿真方法及系统，通过对状态变量的初始值和网络变量的初始值进行预测，能精确获取状态变量的初始值和网络变量的初始值可以有效降低一个时间步长中的迭代次数；在相邻两次迭代的雅可比矩阵解的残差小于第一预设阈值时，不更新雅可比矩阵，减少了更新雅克比矩阵的占用时间，缩短仿真计算时间；本发明实施例通过采用牛顿积分算法求解微分-代数方程组，确保了仿真的精度，在迭代过程中不更新雅克比矩阵，能够缩短机电暂态过程的仿真计算时间，满足实时监控电力系统的需求。
以上所述是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也视为本发明的保护范围。

Claims (10)

1.一种机电暂态过程的仿真方法，其特征在于，包括：
利用电力系统的元件参数和网络拓扑结构，形成描述所述电力系统机电暂态过程的第一微分-代数方程组；
利用所述电力系统初始化时的运行参数，获得初始仿真结果，其中，所述初始仿真结果包括状态变量的初始值和网络变量的初始值；
基于牛顿积分算法，根据所述初始仿真结果计算所述第一微分-代数方程组的雅可比矩阵，并在相邻两次迭代的雅可比矩阵解的残差小于第一预设阈值时，不更新雅可比矩阵，获得所述雅可比矩阵、所述第一微分-代数方程组的第一状态变量和所述第一微分-代数方程组的第一网络变量；
利用所有设备的注入电流，形成描述所述电力系统机电暂态过程的第二微分-代数方程组；
基于牛顿积分算法，根据所述雅可比矩阵、所述第一状态变量和所述第一网络变量对所述第二微分-代数方程组进行迭代求解，在相邻两次迭代的状态变量的残差小于第二预设阈值，且相邻两次迭代的网络变量的残差小于第三预设阈值时，获得所述第二微分-代数方程组的第二状态变量和所述第二微分-代数方程组的第二网络变量；
当仿真时间达到仿真终止时间，输出仿真结果；其中，所述仿真结果包括所述电力系统各个节点的电压、功率、发电机功角曲线和发电机相对摇摆角曲线。
2.如权利要求1所述的机电暂态过程的仿真方法，其特征在于，所述第一微分-代数方程组如公式(1)所示：
Figure FDA0003031259980000011
其中，微分方程表示电力系统元件的动态特性，是电力系统的状态方程；代数方程表示电力系统元件的静态特性，主要是系统的网络方程；x为n个状态向量，V为m个代数向量，微分代数方程组的阶数为n和m之和；Y是电力系统的导纳矩阵。
3.如权利要求2所述的机电暂态过程的仿真方法，其特征在于，所述利用所述电力系统初始化时的运行参数，获得初始仿真结果，其中，所述初始仿真结果包括状态变量的初始值和网络变量的初始值，具体为：
利用所述电力系统初始化时的运行参数，采用简单欧拉法对初始状态变量进行预测，获得状态变量的初始值，则所述状态变量的初始值
Figure FDA0003031259980000021
为：
Figure FDA0003031259980000022
其中，h为仿真步长，xn为第n步的状态变量，Vn为第n步的网络变量；
利用所述电力系统初始化时的运行参数，采用几何预测法对初始网络变量进行预测，获得网络变量的初始值，则所述网络变量的初始值
Figure FDA0003031259980000023
为：
Figure FDA0003031259980000024
其中，Vn-1为第n-1步的网络变量，Vn为第n步的网络变量。
4.如权利要求2所述的机电暂态过程的仿真方法，其特征在于，所述基于牛顿积分算法，根据所述初始仿真结果计算所述第一微分-代数方程组的雅可比矩阵，具体为：
对所述第一微分-代数方程组采用隐式梯形积分法进行差分，得到公式(2)和(3)：
Figure FDA0003031259980000025
YVn+1＝I(xn+1,Vn+1) (3)
其中，Y是电力系统的导纳矩阵，h是步长；
令公式(2)改成为公式(4)：
Figure FDA0003031259980000026
根据所述状态变量的初始值和所述网络变量的初始值，采用牛顿积分算法求解F(xn+1,Vn+1)＝0，计算所述第一微分-代数方程组的雅可比矩阵。
5.如权利要求4所述的机电暂态过程的仿真方法，其特征在于，所述利用所有设备的注入电流，形成描述所述电力系统机电暂态过程的第二微分-代数方程组，具体为：
对公式(2)和公式(3)采用牛顿积分算法求解，得到公式(10)和(14)：
Figure FDA0003031259980000031
Figure FDA0003031259980000032
其中，AG、BG、CG和YG是雅克比矩阵，雅可比矩阵AG、BG、CG和YG是(x，V)的函数；
利用所有设备的注入电流，得到公式(18)：
Figure FDA0003031259980000033
将公式(18)代替公式(14)，形成所述第二微分-代数方程组，所述第二微分-代数方程组如公式(10)和公式(18)所示：
Figure FDA0003031259980000034
Figure FDA0003031259980000035
其中，YGn是电力系统所有设备对应于非凸极效应的常数矩阵。
6.如权利要求1所述的机电暂态过程的仿真方法，其特征在于，所述相邻两次迭代的状态变量的残差具体为：
Figure FDA0003031259980000036
其中，
Figure FDA0003031259980000037
为第k+1次迭代的状态变量，
Figure FDA0003031259980000038
为第k次迭代的状态变量；
所述相邻两次迭代的网络变量的残差具体为：
Figure FDA0003031259980000039
其中，
Figure FDA00030312599800000310
为第k+1次迭代的状态变量，
Figure FDA00030312599800000311
第k次迭代的状态变量。
7.如权利要求1所述的机电暂态过程的仿真方法，其特征在于，所述方法还包括：
当仿真时间未达到仿真终止时间，继续进行下一步长仿真计算。
8.一种机电暂态过程的仿真系统，包括：
第一仿真模型单元：用于利用电力系统的元件参数和网络拓扑结构，形成描述所述电力系统机电暂态过程的第一微分-代数方程组；
初始化单元：用于利用所述电力系统初始化时的运行参数，获得初始仿真结果，其中，所述初始仿真结果包括状态变量的初始值和网络变量的初始值；
雅克比矩阵单元：用于基于牛顿积分算法，根据所述初始仿真结果计算所述第一微分-代数方程组的雅可比矩阵，并在相邻两次迭代的雅可比矩阵解的残差小于第一预设阈值时，不更新雅可比矩阵，获得所述雅可比矩阵、所述第一微分-代数方程组的第一状态变量和所述第一微分-代数方程组的第一网络变量；
第二仿真模型单元：用于利用所有设备的注入电流，形成描述所述电力系统机电暂态过程的第二微分-代数方程组；
求解单元：用于基于牛顿积分算法，根据所述雅可比矩阵、所述第一状态变量和所述第一网络变量对所述第二微分-代数方程组进行迭代求解，在相邻两次迭代的状态变量的残差小于第二预设阈值，且相邻两次迭代的网络变量的残差小于第三预设阈值时，获得所述第二微分-代数方程组的第二状态变量和所述第二微分-代数方程组的第二网络变量；
输出单元：用于当仿真时间达到仿真终止时间，输出仿真结果；其中，所述仿真结果包括所述电力系统各个节点的电压、功率、发电机功角曲线和发电机相对摇摆角曲线。
9.如权利要求8所述的机电暂态过程的仿真系统，其特征在于，所述第一微分-代数方程组如公式(1)所示：
Figure FDA0003031259980000041
其中，微分方程表示电力系统元件的动态特性，是电力系统的状态方程；代数方程表示电力系统元件的静态特性，主要是系统的网络方程；x为n个状态向量，V为m个代数向量，微分代数方程组的阶数为n和m之和；Y是电力系统的导纳矩阵。
10.如权利要求8所述的机电暂态过程的仿真系统，其特征在于，所述相邻两次迭代的状态变量的残差具体为：
Figure FDA0003031259980000051
其中，
Figure FDA0003031259980000052
为第k+1次迭代的状态变量，
Figure FDA0003031259980000053
为第k次迭代的状态变量；
所述相邻两次迭代的网络变量的残差具体为：
Figure FDA0003031259980000054
其中，
Figure FDA0003031259980000055
为第k+1次迭代的状态变量，
Figure FDA0003031259980000056
为第k次迭代的状态变量。
CN202110434097.7A 2021-04-21 2021-04-21 一种机电暂态过程的仿真方法及系统 Active CN113221298B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110434097.7A CN113221298B (zh) 2021-04-21 2021-04-21 一种机电暂态过程的仿真方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110434097.7A CN113221298B (zh) 2021-04-21 2021-04-21 一种机电暂态过程的仿真方法及系统

Publications (2)

Publication Number Publication Date
CN113221298A true CN113221298A (zh) 2021-08-06
CN113221298B CN113221298B (zh) 2023-02-24

Family

ID=77088348

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110434097.7A Active CN113221298B (zh) 2021-04-21 2021-04-21 一种机电暂态过程的仿真方法及系统

Country Status (1)

Country Link
CN (1) CN113221298B (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115719955A (zh) * 2022-11-15 2023-02-28 南方电网数字电网研究院有限公司 电力系统机电暂态超大规模微分代数方程联合求解方法
CN115718986A (zh) * 2022-10-31 2023-02-28 南方电网数字电网研究院有限公司 一种基于分布式内存架构的多核并行时域仿真方法

Citations (14)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101446991A (zh) * 2008-08-15 2009-06-03 中国电力科学研究院 一种电力系统全过程动态仿真的数值积分方法
CN102054095A (zh) * 2010-11-26 2011-05-11 天津大学 一种适用于分布式发电系统的积分方法
CN102270845A (zh) * 2011-08-01 2011-12-07 朱建国 一种电力系统暂态稳定数字仿真计算模型及其算法
CN102609575A (zh) * 2012-01-19 2012-07-25 浙江大学 一种基于隐式数值积分的电力系统暂态稳定仿真方法
CN102664397A (zh) * 2012-03-23 2012-09-12 浙江大学 一种基于隐式精细数值积分的电力系统暂态稳定仿真方法
CN103440374A (zh) * 2013-08-23 2013-12-11 中国南方电网有限责任公司电网技术研究中心 基于状态变量节点分析混合框架的电磁暂态仿真建模方法
CN104732033A (zh) * 2015-03-31 2015-06-24 河南行知专利服务有限公司 一种电力系统机电暂态过程的仿真分析方法
CN106295001A (zh) * 2016-08-10 2017-01-04 华北电力大学 适用于电力系统中长时间尺度的准稳态变步长仿真方法
CN106887844A (zh) * 2017-04-20 2017-06-23 清华大学 一种基于稀疏恢复的快速分解配电网三相状态估计方法
CN108462181A (zh) * 2018-01-31 2018-08-28 天津大学 考虑稀疏性的智能配电网潮流雅可比矩阵鲁棒估计方法
CN109149665A (zh) * 2017-06-28 2019-01-04 中国电力科学研究院 用于柔性直流交流电网联合仿真的多速率仿真方法及系统
CN109698511A (zh) * 2019-01-25 2019-04-30 广西电网有限责任公司电力科学研究院 一种基于伴随灵敏度的暂态稳定预防控制优化算法
US20190243320A1 (en) * 2018-02-05 2019-08-08 Mitsubishi Electric Research Laboratories, Inc. Method and Apparatus for Preconditioned Predictive Control
CN110135031A (zh) * 2019-04-30 2019-08-16 东南大学 基于半隐式龙格库塔法的电力系统暂态稳定计算方法

Patent Citations (14)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101446991A (zh) * 2008-08-15 2009-06-03 中国电力科学研究院 一种电力系统全过程动态仿真的数值积分方法
CN102054095A (zh) * 2010-11-26 2011-05-11 天津大学 一种适用于分布式发电系统的积分方法
CN102270845A (zh) * 2011-08-01 2011-12-07 朱建国 一种电力系统暂态稳定数字仿真计算模型及其算法
CN102609575A (zh) * 2012-01-19 2012-07-25 浙江大学 一种基于隐式数值积分的电力系统暂态稳定仿真方法
CN102664397A (zh) * 2012-03-23 2012-09-12 浙江大学 一种基于隐式精细数值积分的电力系统暂态稳定仿真方法
CN103440374A (zh) * 2013-08-23 2013-12-11 中国南方电网有限责任公司电网技术研究中心 基于状态变量节点分析混合框架的电磁暂态仿真建模方法
CN104732033A (zh) * 2015-03-31 2015-06-24 河南行知专利服务有限公司 一种电力系统机电暂态过程的仿真分析方法
CN106295001A (zh) * 2016-08-10 2017-01-04 华北电力大学 适用于电力系统中长时间尺度的准稳态变步长仿真方法
CN106887844A (zh) * 2017-04-20 2017-06-23 清华大学 一种基于稀疏恢复的快速分解配电网三相状态估计方法
CN109149665A (zh) * 2017-06-28 2019-01-04 中国电力科学研究院 用于柔性直流交流电网联合仿真的多速率仿真方法及系统
CN108462181A (zh) * 2018-01-31 2018-08-28 天津大学 考虑稀疏性的智能配电网潮流雅可比矩阵鲁棒估计方法
US20190243320A1 (en) * 2018-02-05 2019-08-08 Mitsubishi Electric Research Laboratories, Inc. Method and Apparatus for Preconditioned Predictive Control
CN109698511A (zh) * 2019-01-25 2019-04-30 广西电网有限责任公司电力科学研究院 一种基于伴随灵敏度的暂态稳定预防控制优化算法
CN110135031A (zh) * 2019-04-30 2019-08-16 东南大学 基于半隐式龙格库塔法的电力系统暂态稳定计算方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
LIGANG ZHAO,ETC.: "Fast Simulation of MMC-HVDC Systems Using A Generalized Equivalent Circuit Model of MMC", 《2020 IEEE 16TH INTERNATIONAL CONFERENCE ON CONTROL & AUTOMATION (ICCA)》 *

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115718986A (zh) * 2022-10-31 2023-02-28 南方电网数字电网研究院有限公司 一种基于分布式内存架构的多核并行时域仿真方法
CN115718986B (zh) * 2022-10-31 2023-12-12 南方电网数字电网研究院有限公司 一种基于分布式内存架构的多核并行时域仿真方法
CN115719955A (zh) * 2022-11-15 2023-02-28 南方电网数字电网研究院有限公司 电力系统机电暂态超大规模微分代数方程联合求解方法

Also Published As

Publication number Publication date
CN113221298B (zh) 2023-02-24

Similar Documents

Publication Publication Date Title
CN113221298B (zh) 一种机电暂态过程的仿真方法及系统
Sjöberg et al. Initializing Wiener–Hammerstein models based on partitioning of the best linear approximation
Transtrum et al. Measurement-directed reduction of dynamic models in power systems
CN114580647B (zh) 量子系统的模拟方法、计算设备、装置及存储介质
Al’shin et al. Rosenbrock schemes with complex coefficients for stiff and differential algebraic systems
CN107069696A (zh) 一种电力系统状态估计的并行计算方法
EP2864911A1 (en) Method for automated assistance to design nonlinear analog circuit with transient solver
CN111859835A (zh) 一种电路互连网络模型的降阶方法、降阶装置及降阶设备
CN111506873A (zh) 基于扩展卡尔曼滤波算法的非线性系统参数估计方法
Dolin et al. Modified methods of circuit simulation of radio engineering devices in the time domain
CN106856327A (zh) 一种含小阻抗支路电网串联补偿快速分解法潮流计算方法
Anghel et al. Stability analysis of power systems using network decomposition and local gain analysis
Stiasny et al. Solving differential-algebraic equations in power systems dynamics with neural networks and spatial decomposition
CN115455745B (zh) 自适应频点采样的扫频方法、系统及相关设备
Zhang et al. Convergence rate analysis of Gaussian belief propagation for Markov networks
Gonzalez Villarreal et al. Shifting strategy for efficient block‐based non‐linear model predictive control using real‐time iterations
CN108879669B (zh) 基于低阶iigd算法的时滞电力系统特征值分析方法
Pepona et al. Identification of piecewise affine LFR models of interconnected systems
Idema et al. Fast Newton load flow
CN111353635B (zh) 基于粒子群算法的传播矩阵模量的优化拟合方法及系统
Lantsov et al. Algorithms for software tools with harmonic balance method
Voss et al. Model order reduction for nonlinear differential algebraic equations in circuit simulation
US20060100830A1 (en) Moment computations of nonuniform distributed coupled RLC trees with applications to estimating crosstalk noise
Florez et al. Numerical methods coupled with Richardson extrapolation for computation of transient power systems
Ponalagusamy et al. A new fourth order embedded RKAHeM (4, 4) method with error control on multilayer raster cellular neural network

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant