CN110929405B - 一种计及风电机组及燃气轮机组的电-气动态分析方法 - Google Patents
一种计及风电机组及燃气轮机组的电-气动态分析方法 Download PDFInfo
- Publication number
- CN110929405B CN110929405B CN201911189353.XA CN201911189353A CN110929405B CN 110929405 B CN110929405 B CN 110929405B CN 201911189353 A CN201911189353 A CN 201911189353A CN 110929405 B CN110929405 B CN 110929405B
- Authority
- CN
- China
- Prior art keywords
- gas
- network
- power
- node
- flow
- 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.)
- Active
Links
- 238000004458 analytical method Methods 0.000 title claims abstract description 12
- VNWKTOKETHGBQD-UHFFFAOYSA-N methane Chemical compound C VNWKTOKETHGBQD-UHFFFAOYSA-N 0.000 claims abstract description 132
- 239000007789 gas Substances 0.000 claims abstract description 107
- 239000003345 natural gas Substances 0.000 claims abstract description 66
- 238000000034 method Methods 0.000 claims abstract description 36
- 230000003068 static effect Effects 0.000 claims abstract description 25
- 238000004364 calculation method Methods 0.000 claims abstract description 17
- 230000005624 perturbation theories Effects 0.000 claims abstract description 10
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 9
- 239000000446 fuel Substances 0.000 claims description 21
- 230000007423 decrease Effects 0.000 claims description 12
- 238000007665 sagging Methods 0.000 claims description 9
- 230000008878 coupling Effects 0.000 claims description 8
- 238000010168 coupling process Methods 0.000 claims description 8
- 238000005859 coupling reaction Methods 0.000 claims description 8
- 230000008859 change Effects 0.000 claims description 6
- 238000002485 combustion reaction Methods 0.000 claims description 6
- 230000005484 gravity Effects 0.000 claims description 6
- 230000001052 transient effect Effects 0.000 claims description 4
- 238000009529 body temperature measurement Methods 0.000 claims description 3
- 230000006835 compression Effects 0.000 claims description 3
- 238000007906 compression Methods 0.000 claims description 3
- 230000003750 conditioning effect Effects 0.000 claims description 3
- 238000010438 heat treatment Methods 0.000 claims description 3
- 230000001133 acceleration Effects 0.000 claims description 2
- 230000008569 process Effects 0.000 abstract description 3
- 238000004088 simulation Methods 0.000 description 5
- 238000006243 chemical reaction Methods 0.000 description 3
- 230000007547 defect Effects 0.000 description 3
- 230000033228 biological regulation Effects 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000006467 substitution reaction Methods 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000005684 electric field Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000010248 power generation Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000009466 transformation Effects 0.000 description 1
- 238000012795 verification Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q50/00—Information and communication technology [ICT] specially adapted for implementation of business processes of specific business sectors, e.g. utilities or tourism
- G06Q50/06—Energy or water supply
-
- H—ELECTRICITY
- H02—GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
- H02J—CIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
- H02J3/00—Circuit arrangements for ac mains or ac distribution networks
-
- H—ELECTRICITY
- H02—GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
- H02J—CIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
- H02J3/00—Circuit arrangements for ac mains or ac distribution networks
- H02J3/04—Circuit arrangements for ac mains or ac distribution networks for connecting networks of the same frequency but supplied from different sources
- H02J3/06—Controlling transfer of power between connected networks; Controlling sharing of load between connected networks
-
- 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
- Y02E10/00—Energy generation through renewable energy sources
- Y02E10/70—Wind energy
- Y02E10/76—Power conversion electric or electronic aspects
Landscapes
- Engineering & Computer Science (AREA)
- Business, Economics & Management (AREA)
- Power Engineering (AREA)
- Health & Medical Sciences (AREA)
- Economics (AREA)
- Human Resources & Organizations (AREA)
- Tourism & Hospitality (AREA)
- General Health & Medical Sciences (AREA)
- Public Health (AREA)
- Marketing (AREA)
- Primary Health Care (AREA)
- Strategic Management (AREA)
- Water Supply & Treatment (AREA)
- Physics & Mathematics (AREA)
- General Business, Economics & Management (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Supply And Distribution Of Alternating Current (AREA)
- Control Of Eletrric Generators (AREA)
Abstract
本发明提供一种计及风电机组及燃气轮机组的电‑气动态分析方法,涉及综合能源系统技术领域。本发明的主要步骤为:输入电力系统和天然气系统的基本参数,建立电‑气网络静态潮流模型,分别建立网络、风电机组、燃气轮机以及电力系统的动态模型,在动态模型中确定合适的小参量,并基于奇异摄动理论建立系统的多时间尺度算法,采用多时间尺度算法求解网络动态潮流模型。本发明在计算过程中考虑了电‑气互联网络中存在的多时间尺度问题,还可以准确计算电力系统在大扰动后的全网潮流分布;通过使用动态建模的方法,使得所提出的计算方法更为准确。
Description
技术领域
本发明涉及综合能源系统技术领域,尤其涉及一种计及风电机组及燃气轮机组的电-气动态分析方法。
背景技术
为应对能源危机和环境问题,以多种能源耦合为核心的能源互联网应运而生。随着电转气技术的不断发展以及燃气轮机的数量不断提升,电-气网络之间的耦合越来越紧密,这影响着电力系统稳定和经济运行。针对于规模不断扩大的电-气网络,系统潮流计算可以获取全网的运行状态,是进行系统规划、运行控制、能源交易的依据和基准。在文献[1]Qing Zeng,Jiakun Fang.Steady-state analysis of the integrated natural gas andelectric power system with bi-directional energy conversion[J].AppliedEnergy,2016,184:1483-1492.所提出的静态潮流计算方法中,并没有考虑网络的多时间尺度特性,也无法准确获得电力系统在大扰动后(短路、负荷增减以及断线故障)的潮流分布,在一定情况下,运行人员无法通过这种方法准确获取的系统运行状况。而文献[2]MaciejChaczykowski.Transient flow in natural gas pipeline-The effect of pipelinethermal model[J].Applied Mathematical Modelling,2010,34:1052-1067.和文献[3]Ruohan Wang,Yuzheng Xie1,Hengxu Zhang,Changgang Li,Wei Li,VladimirTerzija.Dynamic power flow algorithm considering frequency regulation of windpower generators[J].IET Renewable Power Generation,2017,11(8):1218-1225.仅仅是停留在流体力学领域和电气领域,无法计及系统之间的相互影响,并不适用于能源系统级的潮流计算。此外,风机的大规模并网,也给电力系统稳定性带来一定的挑战。
基于之前研究的不足,本发明提出一种动态潮流计算方法,在计算过程中考虑了电-气互联网络中存在的多时间尺度问题,还可以准确计算电力系统在大扰动后全网的潮流分布;通过使用动态建模的方法,使得所提出的计算方法更为准确。
发明内容
针对现有技术的不足,本发明提供一种计及风电机组及燃气轮机组的电-气动态分析方法,发明方案如下:
一种计及风电机组及燃气轮机组的电-气动态分析方法,包括以下步骤:
步骤1、输入电力系统和天然气系统的基本参数,包括电力系统支路参数、变压器参数、发电机输出功率、负荷参数、节点类型、风速、空气密度、风叶扫过的面积、风速、风能的利用系数、初始减载系数、风速的最大功率点以及下垂系数、天然气网络中的气体状态、气源压力初值、输气管道参数、负荷流量以及节点类型,建立电-气网络静态潮流模型;
所述电-气网络静态潮流模型包括电力系统模型、燃气系统模型以及耦合元件模型,其中所述电力系统静态模型为:
式中Pi和Qi为流入节点i的有功功率和无功功率,Ui为流入节点i的电压,Uj为与流入节点i相连的节点j的电压,Gij为节点i和节点j之间的线路导纳,Bij为节点i和节点j之间的线路电纳,δij为节点i和节点j之间的相角差,m为电力系统网络的节点数;
所述天然气网络的静态模型为:
式中fij为流经节点i和节点j的流量,sign(pi,pj)为管道天然气的流动方向,当pi>pj时,气体是由节点i流向节点j,流量fij为正值,反之流量fij为负值;
公式中管道常数Kij为:
其中,z为气体压缩系数,T为气体温度,Lij为管道的长度,δ为天然气比重,Dij为管道直径,Fij为管道摩擦力;
所述耦合元件模型为:
式中fi为燃气轮机所需流量,GHV为天然气的热值,PGTi为燃气轮机所发出的功率。
步骤2、分别建立天然气网络、风电机组、燃气轮机以及电力系统的动态模型;
所述天然气网络的动态模型为:
p=ρ·c2
式中,ρ表示天然气气体的密度,p表示管网内的气体压强,v表示气体在管道中的流速,x表示管道的位置量,t表示时间,T表示气体温度,R表示气体的常数,Dij表示天然气管道的直径,G为重力加速度,f表示管道内壁的摩擦系数,θ表示管道的倾斜程度,c表示声速;
本实施例中在建模过程中假设天然气在管道中等温传输并且管道和外界没有热量的交换,把方程M=ρ·v·A带入管道的偏微分方程中可以得到简化之后的天然气管道方程:
式中A表示管道的横截面积,M表示流量,∣M∣表示流量的绝对值,当流量的流经方向和假设一致时,绝对值取正号,反之为负号。
为解决上述偏微分方程,本实施例使用Kiuchi隐式方法对方程进行离散,使用的这种近似等效的方法为:
经过上述变换之后的方程为:
其中,公式中M>0,并且参数α为:
经过上述离散之后的偏微分方程中Δx和Δt的取值是任意的,即方程是无条件收敛的。
所述风电机组动态模型为:
其中ρair表示空气密度,A表示风叶所扫过的面积,vair表示风速,vcut-in表示风机的切入风速,vmax表示最大风速,Cp表示对应风速下风能的利用系数,k1表示电网频率为额定值下的初始减载水平,Pop为对应风速的最大功率点,k2为下垂系数,其中△Pw表示为:
其中Tw表示风电机组的转动惯量;△Pwi为风电机组功率改变量;△f表示系统频率波动;k表示调频系数;
在此风速段下,风机转动惯惯量的表达式为:
其中,J表示风电机组转动惯量,SN风电机组表示额定容量,w0表示初始减载水平下的转速,v0表示初始减载水平下对应的风速,wop表示风机在最大功率时输出的转速,vop为对应风速;
所述燃气轮机的动态模型包括速度控制环节:
温度控制环节:
燃料控制环节:
压缩机涡轮环节:
式中,W、Y和X表示速度控制环节系数,w和wref表示发电机转速和参考转速,FD1表示速度控制环节输出的参考燃料值,x1速度控制环节中间变量,x2为热敏电阻中间变量,x3温度传感器中间变量,x4为燃气轮机出口温度中间变量,Vce表示温度控制环节的中间变量,FD2表示温度控制系统所选出的燃料信号,K1燃料控制系统系数,wmin为燃气轮机输入的最小燃料值,Tf表示燃料调节系统时间常数,Wf表示进入燃烧室的燃料值,TCR表示燃烧室延时时间常数,TED和TCD分别表示涡轮环节和温度测量环节时间常数,Wf2表示涡轮输出的高温燃气流量,af1、bf1、af2、bf2和cf2表示涡轮环节中函数f1和f2中的变量,Tm表示输出转矩;
所述电力系统的动态模型是指,当系统出现不平衡功率时,由于转子惯性导致系统频率偏差,表示为:
式中,Tj表示所有发电机组等效为一台发电机时的系统的惯性常数,△f表示系统的频率偏差,Pacc表示系统出现的不平衡功率,K表示系统常数,PLi表示系统负荷量;
当系统出现频率偏差时,各个装有调速器的发电机、燃气轮机和风电机组随之改变输出功率,对于发电机有:
TGi为发电机及风电机组的转动惯量,KGi表示发电机组和风电机组有功-频率系数,△PGi表示发电机组和风电机组增发或者减少的有功功率。
步骤3、在网络、风电机组和燃气轮机的动态模型中确定合适的小参量,并基于奇异摄动理论建立系统的多时间尺度算法;
所述奇异摄动理论表示为:
其中x1和y分别为慢时间尺度变量和快时间尺度变量,ε表示系统的小参量;当系统进入慢时间尺度中,取ε=0,有:
ys表示快时间尺度变量在慢时间尺度下的状态;
当系统进入快时间尺度中,有:
yf表示慢时间尺度变量在快时间尺度下的状态;
取ε=0时,Ω3系统变为:
公式中ξ和T分别表示慢时间尺度步长和快时间尺度步长,显然在BLS系统中,慢时间尺度进入稳态;根据上述奇异摄动理论,在动态模型中确定合适的小参量ε为:
ε=[T+Tw,TGα,Twβ,Y,T3,T4,b,TCD,TCR,TED]T;
其中α和β分别表示发电机组的个数和风电机组的个数。
步骤4、采用多时间尺度算法求解网络动态潮流模型。
步骤4.1:对网络进行初始化,输入电力系统支路参数、变压器参数、发电机输出功率、负荷参数、节点类型、惯性常数、系统有功-频率特性;输入风机参数包括风速、空气密度、风叶扫过的面积、风速、风能的利用系数、初始减载系数、风速的最大功率点和下垂系数;输入天然气网络中的气体状态、气源压力初值、输气管道参数、负荷流量以及节点类型;输入燃气轮机各个环节的参数,并随机设置电力系统和天然气的网络初值Xe0、Xg0以及收敛精度δ;
步骤4.2:对电力系统进行静态潮流计算,并使初值的变化量△Xe小于所设置的精度,记录计算所得结果Xek;然后通过燃气轮机计算其所需的流量值,并计算天然气网络潮流,使天然气网络的变化量△Xg小于所设置的精度,记录计算所得到的各个节点的压强和管道流量Xgk;
步骤4.3:当电力系统出现负荷的增减时,开始计算动态潮流;
步骤4.4:把静态潮流计算所得到的Xek和Xgk作为动态潮流的初始值Xeg0和Xgg0,并选择合适的小参量ε,构造电-气网络奇异摄动系统,首先构造奇异摄动系统Ω1:
选择合适的小参量ε:
ε=[T+Tw,TGi,Twj,Y,T3,T4,b,TCD,TCR,TED]T;
其中α和β分别表示发电机组的个数和风电机组的个数;
取ε=0,慢时间变量进入暂态,快时间变量进入静态,所求得的系统Ω2为:
所得到的BLS系统为:
使慢时间常数和快时间常数t1=0和T1=0,并设置慢时间尺度步长△t和快时间尺度步长△T;
步骤4.5:更新慢时间常数t1,在慢时间尺度步长△t下,由于BLS系统是由一系列微分方程组成的,需要通过改进欧拉法计算BLS系统,计算方法如下:
对于微分方程:
设置x0=0,此时有中间变量:
其中r=0,1,2……,并把此时所求出的f3的值赋给xr+1,有:
xr+1=f3;
Xr+1作为下一次迭代的初始值进行迭代,直到两次计算的值小于设定收敛精度,迭代结束并记录迭代时间r△T1;结束的标志为:
|xr+1-xr|≤δ;
之后判断t1和T1之间的关系,如果t1>T1,进入步骤4.6;若t1<T1,返回步骤4.5继续通过改进欧拉的方法计算,直到满足条件为止;
步骤4.6:把此时燃气轮机的流量作为天然气网络的初值,通过kiuchi隐式方法离散天然气动态潮流方程,通过牛顿法解决网络潮流方程,记录网络动态潮流结果;如果T1+△T>tmax,则继续判断网络是否出现扰动;如果T1+△T<tmax,摄动系统并没有完全解决,需要进入步骤步骤4.5;
步骤4.7:最终输出动态潮流计算结果。
本发明的有益效果:
针对当前电-气互联网络潮流计算方法的不足,本发明在考虑多时间尺度问题的基础上,构建网络动态潮流模型,并提供一种计及风电机组及燃气轮机组的电-气动态分析方法。当电力系统出现负荷的增减、短路或者断线等故障时,本发明所提供的方法可以给系统运行人员提供快速准确的系统状态信息。
附图说明
图1为本发明电-气动态分析方法流程图;
图2为本发明实施例中电力系统和天然气网络系统示意图;
图3为本发明实施例中电力系统和天然气网络系统调频图。
具体实施方式
下面结合附图和具体实施例,对本发明的具体实施方式作进一步详细描述。实施例给出了详细的实施方式和具体的操作过程。但本发明的保护范围不限于下述的实施例。
一种计及风电机组及燃气轮机组的电-气动态分析方法,如图1所示,包括以下步骤:
步骤1、输入电力系统和天然气系统的基本参数,包括电力系统支路参数、变压器参数、发电机输出功率、负荷参数、节点类型、风速、空气密度、风叶扫过的面积、风速、风能的利用系数、初始减载系数、风速的最大功率点以及下垂系数、天然气网络中的气体状态、气源压力初值、输气管道参数、负荷流量以及节点类型,建立电-气网络静态潮流模型;
本实施例中所采用的网络结构为天然气3节点和IEEE9节点,在原先的基础上,把bus2的发电机替换成燃气轮机,在bus3上并入一个风电场,所参考的节点参数选自文献“Mohammad Reza Bank Tavakoli,Behrooz Vahidi,and Wolfgang Gawlik.AnEducational Guide to Extract the Parameters of Heavy Duty Gas Turbines Modelin Dynamic Studies Based on Operational Data[J].IEEE Transactions on PowerSystem,2009,24[3]:1336-1372.”和文献“Qing Zeng,Jiakun Fang.Steady-stateanalysis of the integrated natural gas and electric power system with bi-directional energy conversion[J].Applied Energy,2016,184:1483-1492.”。在电力系统中还需设置的参数为:系统的基准容量为100MVA,基准频率为60Hz,在节点10接入由50台额定容量为2MW风电机组组成的风电场,2MW变速风电机组的参数为:风轮半径为39.8m,最大风能利用系数0.48,最优叶尖速比8.1,切入风速3.0m/s,切出风速22m/s。设空气密度ρ=1.225kg/m3,vmin=8.8m/s、vmax=11.13m/s。设电网频率合格范围为[49.8,50.2]Hz,风电机组的初始减载水平kd0=0.2,下垂系数Rw=0.05。在动态潮流计算中,在计算电力系统中使用的改进欧拉法的步长取为0.02s;在计算天然气网络中设置单位距离为5000m,并且使用的单位步长为3600s。最后设置燃气轮机的基本参数,具体参数参见文献“W.I.Rowen,Simplified mathematical representations of heavy-duty gas turbines[J].Trans.ASME,J.Eng.Power,1983,105(1):865-869.”。
所述电-气网络静态潮流模型包括电力系统模型、燃气系统模型以及耦合元件模型,其中所述电力系统静态模型为:
式中Pi和Qi为流入节点i的有功功率和无功功率,Ui为流入节点i的电压,Uj为与流入节点i相连的节点j的电压,Gij为节点i和节点j之间的线路导纳,Bij为节点i和节点j之间的线路电纳,δij为节点i和节点j之间的相角差,m为电力系统网络的节点数;
所述天然气网络的静态模型为:
式中fij为流经节点i和节点j的流量,sign(pi,pj)为管道天然气的流动方向,当pi>pj时,气体是由节点i流向节点j,流量fij为正值,反之流量fij为负值;
公式中管道常数Kij为:
其中,z为气体压缩系数,T为气体温度,Lij为管道的长度,δ为天然气比重,Dij为管道直径,Fij为管道摩擦力;
所述耦合元件模型为:
式中fi为燃气轮机所需流量,GHV为天然气的热值,PGTi为燃气轮机所发出的功率。
步骤2、分别建立天然气网络、风电机组、燃气轮机以及电力系统的动态模型,如图2所示;
所述天然气网络的动态模型为:
p=ρ·c2
式中,ρ表示天然气气体的密度,p表示管网内的气体压强,v表示气体在管道中的流速,x表示管道的位置量,t表示时间,T表示气体温度,R表示气体的常数,Dij表示天然气管道的直径,G为重力加速度,f表示管道内壁的摩擦系数,θ表示管道的倾斜程度,c表示声速;
所述风电机组动态模型为:
其中ρair表示空气密度,A表示风叶所扫过的面积,vair表示风速,vcut-in表示风机的切入风速,vmax表示最大风速,Cp表示对应风速下风能的利用系数,k1表示电网频率为额定值下的初始减载水平,Pop为对应风速的最大功率点,k2为下垂系数,其中△Pw表示为:
其中Tw表示风电机组的转动惯量;△Pwi为风电机组功率改变量;△f表示系统频率波动;k表示调频系数;
在此风速段下,风机转动惯惯量的表达式为:
其中,J表示风电机组转动惯量,SN风电机组表示额定容量,w0表示初始减载水平下的转速,v0表示初始减载水平下对应的风速,wop表示风机在最大功率时输出的转速,vop为对应风速;
所述燃气轮机的动态模型包括速度控制环节:
温度控制环节:
燃料控制环节:
压缩机涡轮环节:
式中,W、Y和X表示速度控制环节系数,w和wref表示发电机转速和参考转速,FD1表示速度控制环节输出的参考燃料值,x1速度控制环节中间变量,x2为热敏电阻中间变量,x3温度传感器中间变量,x4为燃气轮机出口温度中间变量,Vce表示温度控制环节的中间变量,FD2表示温度控制系统所选出的燃料信号,K1燃料控制系统系数,wmin为燃气轮机输入的最小燃料值,Tf表示燃料调节系统时间常数,Wf表示进入燃烧室的燃料值,TCR表示燃烧室延时时间常数,TED和TCD分别表示涡轮环节和温度测量环节时间常数,Wf2表示涡轮输出的高温燃气流量,af1、bf1、af2、bf2和cf2表示涡轮环节中函数f1和f2中的变量,Tm表示输出转矩;
所述电力系统的动态模型是指,当系统出现不平衡功率时,由于转子惯性导致系统频率偏差,表示为:
式中,Tj表示所有发电机组等效为一台发电机时的系统的惯性常数,△f表示系统的频率偏差,Pacc表示系统出现的不平衡功率,K表示系统常数,PLi表示系统负荷量;
当系统出现频率偏差时,各个装有调速器的发电机、燃气轮机和风电机组随之改变输出功率,对于发电机有:
TGi为发电机及风电机组的转动惯量,KGi表示发电机组和风电机组有功-频率系数,△PGi表示发电机组和风电机组增发或者减少的有功功率。
步骤3、在网络、风电机组和燃气轮机的动态模型中确定合适的小参量,并基于奇异摄动理论建立系统的多时间尺度算法;
所述奇异摄动理论表示为:
其中x1和y分别为慢时间尺度变量和快时间尺度变量,ε表示系统的小参量;当系统进入慢时间尺度中,取ε=0,有:
ys表示快时间尺度变量在慢时间尺度下的状态;
当系统进入快时间尺度中,有:
yf表示慢时间尺度变量在快时间尺度下的状态;
取ε=0时,Ω3系统变为:
公式中ξ和T分别表示慢时间尺度步长和快时间尺度步长,显然在BLS系统中,慢时间尺度进入稳态;根据上述奇异摄动理论,在动态模型中确定合适的小参量ε为:
ε=[T+Tw,TGα,Twβ,Y,T3,T4,b,TCD,TCR,TED]T;
其中α和β分别表示发电机组的个数和风电机组的个数。
本实施例中所选取的小参量ε为:
ε=[275,80,80,100,15,1,15,2.5,0.05,0.4,0.01,0.04,0.2]T
步骤4、采用多时间尺度算法求解网络动态潮流模型。
步骤4.1:对网络进行初始化,输入电力系统支路参数、变压器参数、发电机输出功率、负荷参数、节点类型、惯性常数、系统有功-频率特性;输入风机参数包括风速、空气密度、风叶扫过的面积、风速、风能的利用系数、初始减载系数、风速的最大功率点和下垂系数;输入天然气网络中的气体状态、气源压力初值、输气管道参数、负荷流量以及节点类型;输入燃气轮机各个环节的参数,并随机设置电力系统和天然气的网络初值Xe0、Xg0以及收敛精度δ;
步骤4.2:对电力系统进行静态潮流计算,并使初值的变化量△Xe小于所设置的精度,记录计算所得结果Xek;然后通过燃气轮机计算其所需的流量值,并计算天然气网络潮流,使天然气网络的变化量△Xg小于所设置的精度,记录计算所得到的各个节点的压强和管道流量Xgk;
步骤4.3:当电力系统出现负荷的增减时,开始计算动态潮流;
步骤4.4:把静态潮流计算所得到的Xek和Xgk作为动态潮流的初始值Xeg0和Xgg0,并选择合适的小参量ε,构造电-气网络奇异摄动系统,首先构造奇异摄动系统Ω1:
选择合适的小参量ε:
ε=[T+Tw,TGi,Twj,Y,T3,T4,b,TCD,TCR,TED]T;
其中α和β分别表示发电机组的个数和风电机组的个数;
取ε=0,慢时间变量进入暂态,快时间变量进入静态,所求得的系统Ω2为:
所得到的BLS系统为:
使慢时间常数和快时间常数t1=0和T1=0,并设置慢时间尺度步长△t和快时间尺度步长△T;
步骤4.5:更新慢时间常数t1,在慢时间尺度步长△t下,由于BLS系统是由一系列微分方程组成的,需要通过改进欧拉法计算BLS系统,计算方法如下:
对于微分方程:
设置x0=0,此时有中间变量:
其中r=0,1,2……,并把此时所求出的f3的值赋给xr+1,有:
xr+1=f3;
Xr+1作为下一次迭代的初始值进行迭代,直到两次计算的值小于设定收敛精度,迭代结束并记录迭代时间r△T1;结束的标志为:
|xr+1-xr|≤δ;
之后判断t1和T1之间的关系,如果t1>T1,进入步骤4.6;若t1<T1,返回步骤4.5继续通过改进欧拉的方法计算,直到满足条件为止;
步骤4.6:把此时燃气轮机的流量作为天然气网络的初值,通过kiuchi隐式方法离散天然气动态潮流方程,通过牛顿法解决网络潮流方程,记录网络动态潮流结果;如果T1+△T>tmax,则继续判断网络是否出现扰动;如果T1+△T<tmax,摄动系统并没有完全解决,需要进入步骤步骤4.5;
步骤4.7:最终输出动态潮流计算结果。
本实施例所进行的仿真验证为:所设置的仿真总时长为150s,在50s时位于8号母线的负荷1减少5MW有功功率;在100s时位于6号母线的负荷3增加6MW的有功功率,由于天然气网络设置的仿真时长为24小时,为统一仿真时长,上述电力系统在50s时刻对应天然气网络第20个小时。由图3所示的仿真结果可以看出,当系统负荷减少(增加)时,燃气轮机的输出功率随之降低(增加)。由于各个发电机配置了相应的调速系统,系统频率增加(减少)。当燃气轮机的输出功率减少(增加)时,导致进入燃气轮机的天然气流量减少(增加),此时这种影响会传递至天然气网络,使得天然气网络的潮流出现变化。除此之外,也可以看出电力系统和天然气网络存在巨大的时间差异。
以上实施例仅用以说明本发明的技术方案,而非对其限。尽管参照前述实施例对本发明进行了详细的说明,本领域的普通技术人员应当理解;其依然可以对前述实施例所记述方案进行修改,或对其中部分或全部技术特征进行等同替换;而这些修改或替换皆应在有专利要求书所确定的保护范围内。
Claims (1)
1.一种计及风电机组及燃气轮机组的电-气动态分析方法,其特征在于:包括以下步骤:
步骤1、输入电力系统和天然气系统的基本参数,包括电力系统支路参数、变压器参数、发电机输出功率、负荷参数、节点类型、风速、空气密度、风叶扫过的面积、风能的利用系数、初始减载系数、风速的最大功率点以及下垂系数、天然气网络中的气体状态、气源压力初值、输气管道参数、负荷流量以及节点类型,建立电-气网络静态潮流模型;
步骤1中所述电-气网络静态潮流模型包括电力系统静态模型、天然气网络的静态模型以及耦合元件模型,其中所述电力系统静态模型为:
式中Pi和Qi为流入节点i的有功功率和无功功率,Ui为流入节点i的电压,Uj为与流入节点i相连的节点j的电压,Gij为节点i和节点j之间的线路导纳,Bij为节点i和节点j之间的线路电纳,δij为节点i和节点j之间的相角差,m为电力系统网络的节点数;
所述天然气网络的静态模型为:
式中fij为流经节点i和节点j的流量,sign(pi,pj)为管道天然气的流动方向,当pi>pj时,气体是由节点i流向节点j,流量fij为正值,反之流量fij为负值;
公式中管道常数Kij为:
其中,z为气体压缩系数,T为气体温度,Lij为管道的长度,δ为天然气比重,Dij为管道直径,Fij为管道摩擦力;
所述耦合元件模型为:
式中fi为燃气轮机所需流量,GHV为天然气的热值,PGTi为燃气轮机所发出的功率;
步骤2、分别建立天然气网络、风电机组、燃气轮机以及电力系统的动态模型;
步骤2中所述天然气网络的动态模型为:
p=ρ·c2
式中,ρ表示天然气气体的密度,p表示管网内的气体压强,v表示气体在管道中的流速,x表示管道的位置量,t表示时间,T表示气体温度,R表示气体的常数,Dij表示天然气管道的直径,G为重力加速度,f表示管道内壁的摩擦系数,θ表示管道的倾斜程度,c表示声速;
所述风电机组动态模型为:
其中ρair表示空气密度,A表示风叶所扫过的面积,vair表示风速,vcut-in表示风机的切入风速,vmax表示最大风速,Cp表示对应风速下风能的利用系数,k1表示电网频率为额定值下的初始减载水平,Pop为对应风速的最大功率点,k2为下垂系数,其中△Pw表示为:
其中Tw表示风电机组的转动惯量;△Pw为风电机组功率改变量;△f表示系统频率波动;k表示调频系数;
在此风速段下,风机转动惯量的表达式为:
其中,J表示风电机组转动惯量,SN风电机组表示额定容量,w0表示初始减载水平下的转速,v0表示初始减载水平下对应的风速,wop表示风机在最大功率时输出的转速,vop为对应风速;
所述燃气轮机的动态模型包括速度控制环节:
温度控制环节:
燃料控制环节:
压缩机涡轮环节:
式中,W、Y和X表示速度控制环节系数,w和wref表示发电机转速和参考转速,FD1表示速度控制环节输出的参考燃料值,x1速度控制环节中间变量,x2为热敏电阻中间变量,x3温度传感器中间变量,x4为燃气轮机出口温度中间变量,Vce表示温度控制环节的中间变量,FD2表示温度控制系统所选出的燃料信号,K1燃料控制系统系数,wmin为燃气轮机输入的最小燃料值,Tf表示燃料调节系统时间常数,Wf表示进入燃烧室的燃料值,TCR表示燃烧室延时时间常数,TED和TCD分别表示涡轮环节和温度测量环节时间常数,Wf2表示涡轮输出的高温燃气流量,af1、bf1、af2、bf2和cf2表示涡轮环节中函数f1和f2中的变量,Tm表示输出转矩;
所述电力系统的动态模型是指,当系统出现不平衡功率时,由于转子惯性导致系统频率偏差,表示为:
式中,T+TW表示所有发电机组等效为一台发电机时的系统的惯性常数,△f表示系统的频率偏差,Pacc表示系统出现的不平衡功率,K表示系统常数,PLi表示系统负荷量;
当系统出现频率偏差时,各个装有调速器的发电机、燃气轮机和风电机组随之改变输出功率,对于发电机有:
TGi为发电机及风电机组的转动惯量,KGi表示发电机组和风电机组有功-频率系数,△PGi表示发电机组和风电机组增发或者减少的有功功率;
步骤3、在天然气网络、风电机组和燃气轮机的动态模型中确定合适的小参量,并基于奇异摄动理论建立系统的多时间尺度算法;
步骤3中所述奇异摄动理论表示为:
其中x1和y分别为慢时间尺度变量和快时间尺度变量,ε表示系统的小参量;当系统进入慢时间尺度中,取ε=0,有:
ys表示快时间尺度变量在慢时间尺度下的状态;
当系统进入快时间尺度中,有:
yf表示慢时间尺度变量在快时间尺度下的状态;
取ε=0时,Ω3系统变为:
公式中ξ和τ分别表示慢时间尺度步长和快时间尺度步长,显然在BLS系统中,慢时间尺度进入稳态;根据上述奇异摄动理论,在动态模型中确定合适的小参量ε为:
ε=[T+Tw,TGα,Twβ,Y,T3,T4,b,TCD,TCR,TED]T;
其中α和β分别表示发电机组的个数和风电机组的个数;
步骤4、采用多时间尺度算法求解网络动态潮流模型;
步骤4.1:对网络进行初始化,输入电力系统支路参数、变压器参数、发电机输出功率、负荷参数、节点类型、惯性常数、系统有功-频率特性;输入风机参数包括风速、空气密度、风叶扫过的面积、风速、风能的利用系数、初始减载系数、风速的最大功率点和下垂系数;输入天然气网络中的气体状态、气源压力初值、输气管道参数、负荷流量以及节点类型;输入燃气轮机各个环节的参数,并随机设置电力系统和天然气的网络初值Xe0、Xg0以及收敛精度δ;
步骤4.2:对电力系统进行静态潮流计算,并使初值的变化量△Xe小于所设置的精度,记录计算所得结果Xek;然后通过燃气轮机计算其所需的流量值,并计算天然气网络潮流,使天然气网络的变化量△Xg小于所设置的精度,记录计算所得到的各个节点的压强和管道流量Xgk;
步骤4.3:当电力系统出现负荷的增减时,开始计算动态潮流;
步骤4.4:把静态潮流计算所得到的Xek和Xgk作为动态潮流的初始值Xeg0和Xgg0,并选择合适的小参量ε,构造电-气网络奇异摄动系统,首先构造奇异摄动系统Ω1:
选择合适的小参量ε:
ε=[T+Tw,TGi,Twj,Y,T3,T4,b,TCD,TCR,TED]T;
其中α和β分别表示发电机组的个数和风电机组的个数;
取ε=0,慢时间变量进入暂态,快时间变量进入静态,所求得的系统Ω2为:
所得到的BLS系统为:
使慢时间常数和快时间常数t1=0和τ1=0,并设置慢时间尺度步长△t和快时间尺度步长△τ;
步骤4.5:更新慢时间常数t1,在慢时间尺度步长△t下,由于BLS系统是由一系列微分方程组成的,需要通过改进欧拉法计算BLS系统,计算方法如下:
对于微分方程:
设置x0=0,此时有中间变量:
其中r=0,1,2……,并把此时所求出的f3的值赋给xr+1,有:
xr+1=f3;
Xr+1作为下一次迭代的初始值进行迭代,直到两次计算的值小于设定收敛精度,迭代结束并记录迭代时间r△τ1;结束的标志为:
|xr+1-xr|≤δ;
之后判断t1和τ1之间的关系,如果t1>τ1,进入步骤4.6;若t1<τ1,返回步骤4.5继续通过改进欧拉的方法计算,直到满足条件为止;
步骤4.6:把此时燃气轮机的流量作为天然气网络的初值,通过kiuchi隐式方法离散天然气动态潮流方程,通过牛顿法解决网络潮流方程,记录网络动态潮流结果;如果τ1+△τ>tmax,则继续判断网络是否出现扰动;如果τ1+△τ<tmax,摄动系统并没有完全解决,需要进入步骤步骤4.5;
步骤4.7:最终输出动态潮流计算结果。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201911189353.XA CN110929405B (zh) | 2019-11-28 | 2019-11-28 | 一种计及风电机组及燃气轮机组的电-气动态分析方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201911189353.XA CN110929405B (zh) | 2019-11-28 | 2019-11-28 | 一种计及风电机组及燃气轮机组的电-气动态分析方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN110929405A CN110929405A (zh) | 2020-03-27 |
CN110929405B true CN110929405B (zh) | 2023-08-04 |
Family
ID=69847442
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201911189353.XA Active CN110929405B (zh) | 2019-11-28 | 2019-11-28 | 一种计及风电机组及燃气轮机组的电-气动态分析方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN110929405B (zh) |
Families Citing this family (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN111695269A (zh) * | 2020-06-23 | 2020-09-22 | 中国电力科学研究院有限公司 | 多时间断面的电-气综合能源系统状态估计方法及系统和装置 |
CN112149318B (zh) * | 2020-11-24 | 2021-02-19 | 清华四川能源互联网研究院 | 综合能源系统的协调状态确定方法及相关装置 |
CN114492958B (zh) * | 2022-01-07 | 2024-06-04 | 中国能源建设集团广东省电力设计研究院有限公司 | 一种基于多能流计算的能源管理方法和装置 |
Citations (12)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US5642000A (en) * | 1993-05-03 | 1997-06-24 | Cornell Research Foundation, Inc. | Method for preventing power collapse in electric power systems |
EP1431899A2 (en) * | 2002-12-16 | 2004-06-23 | Rockwell Automation Technologies, Inc. | Decentralized autonomous control for complex fluid distribution systems |
JP2005151746A (ja) * | 2003-11-18 | 2005-06-09 | Hitachi Ltd | 熱電併給型系統制御方法及び熱電併給型系統制御装置 |
AU2011254068A1 (en) * | 2009-02-04 | 2012-01-19 | Applied Hybrid Energy Pty Ltd | Electrical power generating system |
CN105703369A (zh) * | 2016-02-04 | 2016-06-22 | 马瑞 | 一种多能耦合输配网多目标随机模糊动态最优能量流建模与求解方法 |
CN107086576A (zh) * | 2017-06-02 | 2017-08-22 | 武汉理工大学 | 一种分布式潮流控制器多时间尺度数学模型建立方法 |
CN107291990A (zh) * | 2017-05-24 | 2017-10-24 | 河海大学 | 基于电‑气互联综合能源系统暂态模型的能量流仿真方法 |
CN109255550A (zh) * | 2018-09-30 | 2019-01-22 | 东北电力大学 | 一种综合能源系统的n-1静态安全分析方法 |
CN109742763A (zh) * | 2019-02-14 | 2019-05-10 | 国网江西省电力有限公司电力科学研究院 | 一种考虑含多类型能源结构的配电网双向潮流计算方法 |
CN109800968A (zh) * | 2018-12-29 | 2019-05-24 | 重庆大学 | 考虑天然气系统热力过程的电-气互联系统概率能流分析方法 |
CN110011358A (zh) * | 2019-05-05 | 2019-07-12 | 国网安徽省电力有限公司培训中心 | 一种配电网负荷状态调节控制器 |
CN110262236A (zh) * | 2019-06-20 | 2019-09-20 | 合肥工业大学 | 一种电力电子接口并网系统模型降阶的降阶变量选取方法 |
-
2019
- 2019-11-28 CN CN201911189353.XA patent/CN110929405B/zh active Active
Patent Citations (12)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US5642000A (en) * | 1993-05-03 | 1997-06-24 | Cornell Research Foundation, Inc. | Method for preventing power collapse in electric power systems |
EP1431899A2 (en) * | 2002-12-16 | 2004-06-23 | Rockwell Automation Technologies, Inc. | Decentralized autonomous control for complex fluid distribution systems |
JP2005151746A (ja) * | 2003-11-18 | 2005-06-09 | Hitachi Ltd | 熱電併給型系統制御方法及び熱電併給型系統制御装置 |
AU2011254068A1 (en) * | 2009-02-04 | 2012-01-19 | Applied Hybrid Energy Pty Ltd | Electrical power generating system |
CN105703369A (zh) * | 2016-02-04 | 2016-06-22 | 马瑞 | 一种多能耦合输配网多目标随机模糊动态最优能量流建模与求解方法 |
CN107291990A (zh) * | 2017-05-24 | 2017-10-24 | 河海大学 | 基于电‑气互联综合能源系统暂态模型的能量流仿真方法 |
CN107086576A (zh) * | 2017-06-02 | 2017-08-22 | 武汉理工大学 | 一种分布式潮流控制器多时间尺度数学模型建立方法 |
CN109255550A (zh) * | 2018-09-30 | 2019-01-22 | 东北电力大学 | 一种综合能源系统的n-1静态安全分析方法 |
CN109800968A (zh) * | 2018-12-29 | 2019-05-24 | 重庆大学 | 考虑天然气系统热力过程的电-气互联系统概率能流分析方法 |
CN109742763A (zh) * | 2019-02-14 | 2019-05-10 | 国网江西省电力有限公司电力科学研究院 | 一种考虑含多类型能源结构的配电网双向潮流计算方法 |
CN110011358A (zh) * | 2019-05-05 | 2019-07-12 | 国网安徽省电力有限公司培训中心 | 一种配电网负荷状态调节控制器 |
CN110262236A (zh) * | 2019-06-20 | 2019-09-20 | 合肥工业大学 | 一种电力电子接口并网系统模型降阶的降阶变量选取方法 |
Non-Patent Citations (4)
Title |
---|
"Dynamics investigation of flywheel energy storage system in discharge mode for future energy internet";Wang Ziyu;《2017 IEEE Conference on Energy Internet and Energy System Integration (EI2)》;20171128;全文 * |
"Optimized Heat Storage Control Strategy for Combined Operation of Wind Power and Thermoelectric Storage Boiler";Yanfeng Ge;《2019 IEEE Sustainable Power and Energy Conference (iSPEC)》;20191123;全文 * |
"基于MMC的分布式潮流控制器多时间尺度模型研究";金英雷;《中国优秀硕士学位论文全文数据库》;20190715;全文 * |
"生成风电功率时间序列场景的双向优化技术";黎静华;《中国电机工程学报》;20140605;全文 * |
Also Published As
Publication number | Publication date |
---|---|
CN110929405A (zh) | 2020-03-27 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN110929405B (zh) | 一种计及风电机组及燃气轮机组的电-气动态分析方法 | |
Pope et al. | Effects of stator vanes on power coefficients of a zephyr vertical axis wind turbine | |
CN107291990B (zh) | 基于电-气互联综合能源系统暂态模型的能量流仿真方法 | |
Derakhshan et al. | The comparison of incomplete sensitivities and Genetic algorithms applications in 3D radial turbomachinery blade optimization | |
He et al. | New radial turbine dynamic modelling in a low-temperature adiabatic compressed air energy storage system discharging process | |
CN101245765A (zh) | 变桨距风力发电系统的逆系统鲁棒控制方法 | |
Brasil Junior et al. | On the design of propeller hydrokinetic turbines: the effect of the number of blades | |
CN109886523B (zh) | 一种综合能源网动态模型多速率计算方法 | |
CN112966394B (zh) | 一种水力耦合条件水电机群动态特性的仿真方法及系统 | |
Zhang et al. | Control optimisation for pumped storage unit in micro‐grid with wind power penetration using improved grey wolf optimiser | |
CN110647040B (zh) | 一种综合能源系统的安全控制方法及装置 | |
Guo et al. | Dynamic characteristics of a hydro‐turbine governing system considering draft tube pressure pulsation | |
Okulov et al. | Power properties of two interacting wind turbine rotors | |
Huang et al. | A multi-rate dynamic energy flow analysis method for integrated electricity-gas-heat system with different time-scale | |
Wang et al. | Non‐linear modelling and stability analysis of the PTGS at pump mode | |
Chu et al. | Comparative study of the performances of a bio-inspired flexible-bladed wind turbine and a rigid-bladed wind turbine in centimeter-scale | |
Ajirlo et al. | Development of a wind turbine simulator to design and test micro HAWTs | |
Jurado et al. | Use of ARX algorithms for modelling micro-turbines on the distribution feeder | |
Raach et al. | ℋ∞ controller design for closed-loop wake redirection | |
Yang | Towards the development of a wake meandering model based on neural networks | |
Mendes et al. | A computational fluid dynamics investigation on the axial induction factor of a small horizontal axis wind turbine | |
CN105652692A (zh) | 基于热发电的电厂仪控系统的半实物仿真平台及控制方法 | |
Kazda et al. | Mitigating adverse wake effects in a wind farm using non-optimum operational conditions | |
CN114285037B (zh) | 一种区域电-气综合能源系统控制参数稳定域确定方法 | |
Pachauri et al. | Pitch angle controlling of wind turbine system using proportional-integral/fuzzy logic controller |
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 |