CN113991647A - 一种面向频率响应容量规划的电力系统随机生产模拟方法 - Google Patents

Abstract

一种面向频率响应容量规划的电力系统随机生产模拟方法，结合电力系统的多资源结构组成以及快速频率响应需求，包括：1)模拟系统各机组的运行工况，在风力发电系统中计算风电实际出力状况。在常规发电机组中，计算出力；2)从低阶模型出发确定系统参数与频率动态变化过程和模型聚合之间的联系。基于参数聚合思想聚合含有多机组的风电、火电、水电、燃气轮机组及储能设备的频率响应模型，降低模型阶数。结合多资源机组的出力情况计算系统的失负荷带入聚合模型中，求解频率动态。本发明能够考虑频率动态变化过程，与传统方法相比能够更加精确，能够满足面向快速频率响应的备用容量规划需求，为以低频减载为指标的备用容量规划提供了新思路。

Description

一种面向频率响应容量规划的电力系统随机生产模拟方法

技术领域

本发明属于电网频率稳定控制领域，涉及一种针对于频率响应容量规划的随机生产模拟方法。

背景技术

频率响应能力是电力系统频率控制的第一道防线，一次调频过程是阻止大功率缺额下电力系统低频减载、维持系统稳定的关键阶段。在以火电为主体的传统电力系统，一次调频能力强，快速频率响应备用容量的需求小，其容量规划主要以供需功率偏差为基准，即以热备用、冷备用及负荷备用等，较为粗犷。然而，随着可能源机组的大规模并网，电力系统控制高度电力电子化以及特高压交直流混联电网逐渐成型，我国电力系统的常规机组容量占比下降，电网转动惯量减小、一次调频能力减弱，面临着的大功率缺失故障风险增加，系统频率安全稳定形势日益严峻。因此，在新形势下，为了维持系统的频率安全稳定，电力系统对快速频率响应的备用容量需求增大，有必要将频率响应容量规划进行单独划分，以解决日益严峻的频率稳定问题。

发明内容

针对现有技术存在的上述问题，本发明提供一种针对于频率响应容量规划的随机生产模拟方法。

本发明采用的技术方案为：

一种面向频率响应容量规划的电力系统随机生产模拟方法，所述模拟方法包括以下步骤：

步骤1：建立电力系统状态概率模型，模拟系统运行状况。本发明将其分为两个步骤：风力发电机组的状态概率模型和常规发电机组的状态概率模型，具体步骤如下：

步骤1.1：根据马尔可夫链蒙特卡洛原理，建立风力发电机组的状态概率模型包括风电机组运行状态转移与降额状态相互转换和运行状态与停运状态的相互转换，如图1所示。λ₁为风机系统由运行状态转移为停运状态的转移变化率，λ₂为风机系统由运行状态转移为降额状态的转移变化率，μ₁为风机系统由停运状态转移为运行状态的修复率，μ₂为风机系统由降额状态转移为运行状态的修复率。根据风机的三状态模型以及马尔可夫理论，建立状态转移矩阵T：

基于马尔可夫过程的逼近原理，定义P为状态概率矩阵，I为单位矩阵且P＝[P₁,P₂,P₃]，P₁为风机处于运行状态的概率，P₂为风机处于降额状态的概率，P₃为风机处于停运状态的概率，得出：

即：

对上式机型转置运算，并增加全概率条件，将矩阵式中第一组方程用全概率的条件来替换，得出：

因此，有上述矩阵方程，可以得出风机各状态的概率：

显然，在已知风机各状态概率的情况下，只需要利用蒙特卡洛方法对系统进行随机抽样，确定风力发电系统的运行状态。

系统处在运行状态尚未进入降额状态的持续时间为τ₁₂：

系统处在运行状态尚未进入停运状态的持续时间为τ₁₃：

其中，MTTF为失效前平均持续运行时间，γ₁₂和γ₁₃为区间[0,1]内的均与分布随机数。

系统由降额状态或停运状态转移为正常运行状态的修复时间分别为τ₂₁和τ₃₁：

其中，MTTR为系统的平均修复时间，γ₂₁和γ₃₁为区间[0,1]内的均与分布随机数。

步骤1.2：同样，根据该原理，建立常规机组的状态概率模型，主要为运行转态和停运状态之间的转换，如图2所示。定义FOR为强迫停运率，即当蒙特卡洛随机抽样结果小于该值，系统进入停运状态；当其结果大于该值时，系统处于正常运行状态，其求解方法如下：

其系统在每个状态下的运行时间和风力发电机组的算法相似，正常运行状态下和停运状态下的运行时间分别为τ₁和τ₂，即：

τ₁＝-MTTF lnγ₁ (11)

τ₂＝-MTTRlnγ₂ (12)

其中，γ₁和γ₂为区间[0,1]之间均匀分布的随机数。

步骤1.3:根据以上的模型理论，本发明对其进行编程实现，发现蒙特卡洛法在进行随机抽样的过程中，会出现连续停运的情况。因此，本发明引入故障过滤的算法，程序流程图如图3所示。具体为：

1)第i台机组进入故障过滤检测，检测其在整个周期内，第j个时刻是否有连续故障趋势。若没有，则跳出故障过滤程序。

2)若具有连续故障趋势，搜寻从该时刻开始共有多少个故障时刻。若超出规定时间(整个周期)，结束搜索计算，将该时刻之后，规定时刻之前的状态赋值0；若未超出规定时间，将该时刻之后的所有故障连续时刻状态赋值0。

3)继续检测故障时刻，直到每台机组的所有时刻的故障状态全部检测完毕，结束故障过滤。

步骤2：通过状态概率模型，得出各类型机组的生产状态，现结合其工况建立电力系统输出功率模型，模拟系统的失负荷量。其主要分为3个步骤：建立风力发电输出功率模型、建立常规机组输出功率模型以及建立失负荷量模型。

步骤2.1：建立风力发电输出功率模型。风力发电机组的出力主要与所处风电场的风速密切相关，具有一定程度的线性关系。首先需要对风电场的风速进行模拟。基于ARMA模型的基本思想对未来的风速进行预测，进行风电场的风速模拟。

基于风速的历史数据，对其进行标准化，得出风速平均值μ和方差σ，t时刻的预测风速SW_t的表达式如下：

SW_t＝μ+σy_t (13)

其中，y_t为预测的时间序列值。

因为ARMA模型不仅与历史数据有关，而且与历史的系统扰动密切相关，所以假设系统扰动的白噪声序列{ε_t}的均值为0，方差为

则时间序列y_t的为：

其中，

为自回归系数，θ_i(j＝1,2,3,…,q)为滑动平均系数，p和q分别ARMA模型的自回归阶数和滑动平均阶数。

对于ARMA模型的求解，需要对其进行定阶。采用AIC准则对ARMA模型进行模型的阶数确定。定义AIC准则函数如下：

AIC(p,q)＝logσ²(p,q)+2(p+q+1)/N (15)

式中，σ²(p,q)为预测时间序列的残差方差，N为样本的个数。

在AIC准则的模型中，根据经济学和数学的相关理论，当p和q达到某一数值时，使得AIC的值最小。此时，使得AIC的值达到最小的p和q为ARMA模型的阶数。

在风力发电机组中，其风机出力与风电场的风速具体的线性关系如下图14所示。在图中，v_ci是时风力机组启动时的切入风速，v_r是风力机组正常运行时的额定风速，v_co是风力机组运行时的切出风速，设置目的是为保证风机安全使风机停止工作；P_r是额定风速下风力发电机组正常运行时的额定功率。当风速超过切入风速v_ci且未到达额定风速v_r时，风机的输出功率与风速是非线性关系；当风速超过额定风速v_r但未到达切出风速v_co时，风机保持额定功率一直不变；当风速超过切出风速v_co时，风机停止工作。将上述过程建立为如下的分段函数：

其中，P_w为t时刻风电机组的输出功率，v_t为t时刻风机的预测风速，α、β、γ为风电机组输出功率的待定系数，具体求解算式如下：

步骤2.2：建立常规机组的输出功率模型。考虑的常规发电机组共有三种类型，包括火力发电机组、水力发电机组和燃气轮机发电机组，本发明主要以算例IEEE RTS-79中各机组的系统容量为基准，计算其输出功率。

步骤2.3在电力系统的总输出功率的部分，建立电力系统的出力主要由常规机组和风电机组提供。因此，需要分为常规和风电机组两个部分，将系统的输出功率与系统状态相结合，分别进行计算，即：

其中，P为系统的总输出功率，P_i为各类型机组的输出功率，α_i为各类型机组的系统状态。

因此，电力系统的失负荷量ΔP_d为：

ΔP_d＝P-P_L (19)

其中，P_L为算例中的电力系统负荷。

步骤3：将上述模拟计算出的系统失负荷量带入到系统频率响应模型中，求解系统的频率动态。因此，需建立各资源下系统的频率响应模型。本发明从低阶频率响应模型出发，进行模型的建立与分析。下面将逐个分析建立火电、水电、燃气轮机组以及储能设备的低阶频率响应模型。

步骤3.1：建立火电机组的低阶频率响应模型。火电机组的结构图如下图5所示，机组功率增量高压缸汽轮机惯性环节和机械功率增益环节，与负荷侧扰动一同输入到发电机惯性环节中，输出频差，最后频差经调速器负反馈环节，反馈到输入端。其中,P_e为发电机的负荷侧功率、P_m为汽轮机的机械功率、P_SP为机组功率增量，P_α为系统加速功率。T_R为再热时间常数、H为系统惯性常数、1/R_T负反馈增益(R_T为调差系数)、F_H为高压汽轮机功率比例系数、K_T为机械功率增益、D为阻尼系数、T₁为调速器时间常数、T₂为汽容时间常数。

步骤3.2：建立水电机组的低阶频率响应模型。水电机组的结构图如下图6所示，负荷侧扰动，经水电机组的发电机转子惯性环节，输出频差，最后频差经调速器负反馈环节至输入端。其中，H为机组的惯性时间常数，D为机组的阻尼常数，T_y为接力器的响应时间常数，T_w为水流的惯性时间常数，T_RH为水轮机调速器暂态下垂时间常数，T_RS为水轮机调速器复位时间，1/R_h为水轮机的调差系数。

步骤3.3：建立燃气轮机组的低阶频率响应模型。燃气轮机组的结构图如下图7所示，负荷侧扰动，经燃气轮机机组的发电机转子惯性环节，输出频差，最后频差经调速器负反馈环节至输入端。其中，R_G为系统的调差系数，K_G为功率容量占比，c_g和b_g为燃油系统中阀门定位器的传递函数的系数，X_g是调速器部分的超前时间常数，Y_g是调速器部分的滞后时间常数，T_CR是燃烧反应延迟时间常数，T_F是燃油系统时间常数，T_CD是压缩排量时间常数。

步骤3.4：建立储能设备的低阶频率响应模型。储能设备的结构图如图8所示，其中，T_ESD为储能环节的储能时间常数。

步骤4：基于各资源低阶频率响应模型的分析后，为了提升计算速度，具体为整个系统扰动作为输入量，经聚合发电机惯量，输出为频差频差经各储能调节环节反馈至输入端。本发明建立了多资源、多机系统频率响应聚合模型。建模的核心思想是：先聚合多机，再聚合多资源。

步骤4.1：建立多机聚合的系统频率响应模型，火电机组的多机模型如图9所示，具体为整个系统扰动作为输入量，经聚合发电机惯量，输出为频差，频差经各机组调速器反馈至输入端。在多机系统频率响应模型中，每台机组的调差系数R_ti和机械功率增益K_ti与发电机的视在功率S_i密切相关，其中K_ti表示机组在整个电力系统发电机组的功率占比。因此，为了聚合模型，做了如下定义为：

其中,κ_ti为火电机组的等效增益参数。为了简化聚合模型，对调速器机组做归一化处理，定义每个火电机组支路的标准化增益λ_ti为：

F_H、T_R、T_G、T_C、T_R为聚合后频率响应模型的综合等效参数。在聚合计算中，需要找到火电参数最优解X＝{T_G,T_C,F_H,T_R},即聚合的频率响应模型参数可以等效为每一个机组参数的加权平均值。

水电机组的多机频率响应模型如图6所示，具体为整个系统扰动作为输入量，经聚合发电机惯量，输出为频差，频差经各机组调速器反馈至输入端。定义聚合模型的调差系数为：

其中,κ_hi为水电机组的等效增益参数。定义每个水电机组支路的标准化增益λ_hi为：

令燃气轮机参数最优解Y＝{T_w,T_y,T_RS,T_RH}，其系统参数聚合为：

燃气轮机组的多机频率响应模型如图7所示，具体为整个系统扰动作为输入量，经聚合发电机惯量，输出为频差，频差经各机组调速器反馈至输入端。定义聚合模型调差系数为：

其中,κ_gi为水电机组的等效增益参数。定义每个水电机组支路的标准化增益λ_gi为：

令水电参数最优解Z＝{c_g,b_g,X_g,Y_g,T_CR,T_F,T_CD}，其系统参数聚合为：

步骤4.2：将上述各资源模型多机分别聚合为单机模型后，进行多资源模型聚合。根据各资源模型的聚合理论，可以发现机械功率增益K与系统的容量占比密切相关，即：

其中，S_i为各资源机组的系统容量，S_sys为整个多资源机组的系统容量。

在发电机及负荷部分，以荷侧功率为基准，计算每种类型机组的容量占比为：

其中，H_i为各机组的发电机组的惯性时间常数，H为聚合后系统的发电机组的惯性时间常数。

在含有多资源、多机组电力系统频率响应聚合模型中，由于风电机组不具有实际的转子惯性，本发明不考虑其调频特性。依据“先聚合多机，再聚合多资源”的核心思想，构建多资源、多机系统频率响应聚合模型如图12所示，具体为：整个系统扰动作为输入量，经聚合发电机惯量，输出为频差，频差经各资源机组的聚合调速器反馈至输入端。

步骤5：为了提高计算速率，现对多资源、多机系统频率响应聚合模型进行开环解析。将上述含多种资源的电力系统频率响应聚合模型，进行了开环简化，并将各类型机组的调速器系统均拟合为一阶惯性环节。此外，本发明还忽略负荷阻尼D对频率的影响，其输出量为每台机组的有功出力，模型如图13所示，具体为：以系统负荷扰动为输入量，经发电机惯性环节输出为频差，将开环后频差信号经各资源机组的聚合调速器反馈至输入端。

根据开环模型，每台机组输出功率的计算如式所示。

式中，K_i表示各机组拟合后的功率—频率静特性系数；T_i表示各机组响应时间常数，ΔP_Gi为每台机组的调速器增发功率，Δω为系统的频差。在实际工程应用中，可认为在秒级尺度内系统频率初始衰减的斜率为常数m_Δω如图10所示，即：在扰动时间后，系统频差与时间成线性关系，解析式如下。

由此，对ΔP_G进行拉普拉斯逆变换可得其时域表达式分别为：

根据简化后的开环模型和叠加定理，Δω的时域表达式为：

令调速器拟合参数

积分后为：

当系统频率达到最低点时，有

此时时间为t_nadir，即频率最低点到达时间：

但上述式子中出现的C_i、t_nadir均为未知量，即共有N+1个未知量，故须建立如下的N+1次方程组进行求解。

最大频差Δω_max为：

频率最低点f_min为：

f_min＝f₀-f_B·Δω_max (39)

式中，f₀为系统扰动前稳态频率，f_B为系统基础频率。

步骤6：将上述模拟计算出的系统失负荷量带入到系统频率响应模型中，求解系统的频率动态，其主要动态指标为：频率最大偏差、频率最低点到达时间以及频率的初始变化率。

步骤7：通过MATLAB方式实现上述过程，其仿真算例及其结果表明：本发明提出随机生产模拟方法考虑频率动态变化过程，与传统方法相比能够更加精确，能够满足面向快速频率响应的备用容量规划需求，为以低频减载为指标的备用容量规划提供了新思路。

可选择地，含有多资源、多机组的系统频率响应聚合模型，在模型中以单机的形式存在，主要包括：风电、火电、水电、燃气轮机以及储能设备，覆盖面较广，能有效反应现代电力系统现状。此外，聚合模型对于频率动态指标的求解，可以降低模型阶数，有效提升运算速度。

可选择地，引入了大规模的风力发电，并且在传统计算方法上考虑了电力系统的工况，并加入故障过滤的优化算法，更好地满足了实际的电力生产需求。随后，结合相关的风速预测计算方法和算例，模拟系统的出力和失负荷量，具有很强的实际意义。

可选择地，在频率响应容量规划过程中，需考虑频率动态变化过程，来面对在大功率缺额下的快速频率响应需求。

本发明的有益效果是：以最大频率偏差为主要动态指标，与准稳态频率相比更能反映频率的动态变化过程。在面向频率响应过程制定备用容量规划时，考虑该动态指标对频率恶化实行有效拦截具有显著优势。此外，本发明提出将多资源、多机组的频率响应模型聚合并解析，直接对频率最低点求解，降低模型阶数，简化运算过程，极大地提高了计算效率，并且与时序仿真结果相比，其计算精度也在可接受裕度范围。

附图说明

图1为本发明提供的一种风力发电机组的状态转移模型；

图2为本发明提供的一种常规发电机组的状态转移模型；

图3为本发明提供的一种系统故障过滤模型；

图4为本发明提供的一种风力发电机组的输出功率与风速的模型；

图5为本发明提供的一种火电机组的低阶系统频率响应模型；

图6为本发明提供的一种水电机组的低阶系统频率响应模型；

图7为本发明提供的一种燃气机组的低阶系统频率响应模型；

图8为本发明提供的一种储能设备的低阶系统频率响应模型；

图9为本发明提供的一种火电机组的多机系统频率响应模型；

图10为本发明提供的一种水电机组的多机系统频率响应模型；

图11为本发明提供的一种燃气轮机机组的多机系统频率响应模型；

图12为本发明提供的一种含多资源的系统频率响应聚合模型；

图13为本发明提供的一种含多资源的系统频率响应聚合开环模型；

图14为本发明提供的一种系统频率变化初始阶段模型示意图；

图15为本发明提供的火电机组的多机频率响应模型与聚合模型仿真结果；

图16为本发明提供的水电机组的多机频率响应模型与聚合模型仿真结果；

图17为本发明提供的燃气机组的多机频率响应模型与聚合模型仿真结果；

图18为本发明提供的燃气轮机组的惯性环节拟合结果；

图19为本发明提供的一种一年内电力系统失负荷量曲线；

图20为本发明提供的一种准稳态频率、时序仿真法以及解析法的对比。

具体实施方式

下面结合附图和具体的实施对本发明作进一步的说明，但是不作为本发明的限定。

面向频率响应的容量规划，建立了一种考虑频率动态的随机生产模拟方法。结合电力系统的多资源结构组成以及快速频率响应需求，主要将其分为两部分：多资源、多机组系统频率响应模型和电力系统的失负荷量求解。第一部分，本文依据马尔可夫链蒙特卡洛的方法，模拟系统各机组的运行工况。在风力发电系统中，为了贴合实际运行状况，本发明利用历史数据和ARMA模型预测了风电场的风速，结合风电出力函数，计算风电实际出力状况。在常规发电机组中，主要引入改进后的IEEE RTS-79算例进行计算出力。第二部分，首先从低阶模型出发，确定系统参数与频率动态变化过程和模型聚合之间的联系。基于参数聚合的思想，聚合含有多机组的风电、火电、水电、燃气轮机组以及储能设备的频率响应模型，降低模型阶数，有效提升计算速度，便于求解频率动态。最后，结合多资源机组的出力情况计算系统的失负荷带入聚合模型中，求解频率动态与传统方法形成鲜明对比。仿真算例及其结果表明：本文提出随机生产模拟方法考虑频率动态变化过程，与传统方法相比能够更加精确，能够满足面向快速频率响应的备用容量规划需求，为以低频减载为指标的备用容量规划提供了新思路。所述方法设计具体包括以下步骤：

步骤1：建立电力系统状态概率模型，模拟系统运行状况。本发明将其分为两个步骤：风力发电机组的状态概率模型和常规发电机组的状态概率模型，具体步骤如下：

步骤1.1：根据马尔可夫链蒙特卡洛原理，建立风力发电机组的状态概率模型，包括风电机组运行状态转移与降额状态相互转换和运行状态与停运状态的相互转换，如图1所示。λ₁为风机系统由运行状态转移为停运状态的转移变化率，λ₂为风机系统由运行状态转移为降额状态的转移变化率，μ₁为风机系统由停运状态转移为运行状态的修复率，μ₂为风机系统由降额状态转移为运行状态的修复率。根据风机的三状态模型以及马尔可夫理论，建立状态转移矩阵T如公式(1)所示。

基于马尔可夫过程的逼近原理，定义P为状态概率矩阵，I为单位矩阵且P＝[P₁,P₂,P₃]，P₁为风机处于运行状态的概率，P₂为风机处于降额状态的概率，P₃为风机处于停运状态的概率，得出：

即：

对上式机型转置运算，并增加全概率条件，将矩阵式中第一组方程用全概率的条件来替换，得出：

因此，有上述矩阵方程，可以得出风机各状态的概率：

显然，在已知风机各状态概率的情况下，只需要利用蒙特卡洛方法对系统进行随机抽样，确定风力发电系统的运行状态。

系统处在运行状态尚未进入降额状态的持续时间为公式(6)所示的τ₁₂：

系统处在运行状态尚未进入停运状态的持续时间为公式(7)所示的τ₁₃：

系统由降额状态或停运状态转移为正常运行状态的修复时间分别为公式(8)所示的τ₂₁和为公式(9)所示的τ₃₁。

步骤1.2：同样，根据该原理，建立常规机组的状态概率模型，主要为运行转态和停运状态之间的转换，如图2所示。定义FOR为强迫停运率，即当蒙特卡洛随机抽样结果小于该值，系统进入停运状态；当其结果大于该值时，系统处于正常运行状态，其求解方法如下：

其系统在每个状态下的运行时间和风力发电机组的算法相似，正常运行状态下和停运状态下的运行时间分别为τ₁和τ₂，即：

τ₁＝-MTTF lnγ₁ (11)

τ₂＝-MTTRlnγ₂ (12)

其中，γ₁和γ₂为区间[0,1]之间均匀分布的随机数。

步骤1.3:根据以上的模型理论，本发明对其进行编程实现，发现蒙特卡洛法在进行随机抽样的过程中，会出现连续停运的情况。因此，本发明引入故障过滤的算法，程序流程图如图3所示。具体为：

1)第i台机组进入故障过滤检测，检测其在整个周期内，第j个时刻是否有连续故障趋势。若没有，则跳出故障过滤程序。

2)若具有连续故障趋势，搜寻从该时刻开始共有多少个故障时刻。若超出规定时间(整个周期)，结束搜索计算，将该时刻之后，规定时刻之前的状态赋值0；若未超出规定时间，将该时刻之后的所有故障连续时刻状态赋值0。

3)继续检测故障时刻，直到每台机组的所有时刻的故障状态全部检测完毕，结束故障过滤。

步骤2：通过状态概率模型，得出各类型机组的生产状态，现结合其工况建立电力系统输出功率模型，模拟系统的失负荷量。其主要分为3个步骤：建立风力发电输出功率模型、建立常规机组输出功率模型以及建立失负荷量模型。

步骤2.1：建立风力发电输出功率模型。风力发电机组的出力主要与所处风电场的风速密切相关，具有一定程度的线性关系。首先需要对风电场的风速进行模拟。基于ARMA模型的基本思想对未来的风速进行预测，进行风电场的风速模拟。

基于风速的历史数据，对其进行标准化，得出风速平均值μ和方差σ，t时刻的预测风速SW_t的表达式如下：

SW_t＝μ+σy_t (13)

其中，y_t为预测的时间序列值。

因为ARMA模型不仅与历史数据有关，而且与历史的系统扰动密切相关，所以假设系统扰动的白噪声序列{ε_t}的均值为0，方差为

则时间序列y_t的为：

其中，

为自回归系数，θ_i(j＝1,2,3,…,q)为滑动平均系数，p和q分别ARMA模型的自回归阶数和滑动平均阶数。

对于ARMA模型的求解，需要对其进行定阶。采用AIC准则对ARMA模型进行模型的阶数确定。定义AIC准则函数如下：

AIC(p,q)＝logσ²(p,q)+2(p+q+1)/N (15)

式中，σ²(p,q)为预测时间序列的残差方差，N为样本的个数。

在AIC准则的模型中，根据经济学和数学的相关理论，当p和q达到某一数值时，使得AIC的值最小。此时，使得AIC的值达到最小的p和q为ARMA模型的阶数。

在风力发电机组中，其风机出力与风电场的风速具体的线性关系如下图14所示。在图中，v_ci是时风力机组启动时的切入风速，v_r是风力机组正常运行时的额定风速，v_co是风力机组运行时的切出风速，设置目的是为保证风机安全使风机停止工作；P_r是额定风速下风力发电机组正常运行时的额定功率。当风速超过切入风速v_ci且未到达额定风速v_r时，风机的输出功率与风速是非线性关系；当风速超过额定风速v_r但未到达切出风速v_co时，风机保持额定功率一直不变；当风速超过切出风速v_co时，风机停止工作。将上述过程建立为如下的分段函数：

其中，P_w为t时刻风电机组的输出功率，v_t为t时刻风机的预测风速，α、β、γ为风电机组输出功率的待定系数，具体求解算式如下：

步骤2.2：建立常规机组的输出功率模型。考虑的常规发电机组共有三种类型，包括火力发电机组、水力发电机组和燃气轮机发电机组，本发明主要以算例IEEE RTS-79中各机组的系统容量为基准，计算其输出功率。

步骤2.3在电力系统的总输出功率的部分，建立电力系统的出力主要由常规机组和风电机组提供。因此，需要分为常规和风电机组两个部分，将系统的输出功率与系统状态相结合，分别进行计算，即：

其中，P为系统的总输出功率，P_i为各类型机组的输出功率，α_i为各类型机组的系统状态。

因此，电力系统的失负荷量ΔP_d为：ΔP_d＝P-P_L (19)

其中，P_L为算例中的电力系统负荷。

步骤3：将上述模拟计算出的系统失负荷量带入到系统频率响应模型中，求解系统的频率动态。因此，需建立各资源下系统的频率响应模型。本发明从低阶频率响应模型出发，进行模型的建立与分析。下面将逐个分析建立火电、水电、燃气轮机组以及储能设备的低阶频率响应模型。

步骤3.1：建立火电机组的低阶频率响应模型。火电机组的结构图如下图5所示，机组功率增量高压缸汽轮机惯性环节和机械功率增益环节，与负荷侧扰动一同输入到发电机惯性环节中，输出频差，最后频差经调速器负反馈环节，反馈到输入端。其中,P_e为发电机的负荷侧功率、P_m为汽轮机的机械功率、P_SP为机组功率增量，P_α为系统加速功率。T_R为再热时间常数、H为系统惯性常数、1/R_T负反馈增益(R_T为调差系数)、F_H为高压汽轮机功率比例系数、K_T为机械功率增益、D为阻尼系数、T₁为调速器时间常数、T₂为汽容时间常数。

步骤3.2：建立水电机组的低阶频率响应模型。水电机组的结构图如下图6所示，负荷侧扰动，经水电机组的发电机转子惯性环节，输出频差，最后频差经调速器负反馈环节至输入端。其中，其中，H为机组的惯性时间常数，D为机组的阻尼常数，T_y为接力器的响应时间常数，T_w为水流的惯性时间常数，T_RH为水轮机调速器暂态下垂时间常数，T_RS为水轮机调速器复位时间，1/R_h为水轮机的调差系数。

步骤3.3：建立燃气轮机组的低阶频率响应模型。燃气轮机组的结构图如下图7所示，负荷侧扰动，经燃气轮机机组的发电机转子惯性环节，输出频差，最后频差经调速器负反馈环节至输入端。其中，R_G为系统的调差系数，K_G为功率容量占比，c_g和b_g为燃油系统中阀门定位器的传递函数的系数，X_g是调速器部分的超前时间常数，Y_g是调速器部分的滞后时间常数，T_CR是燃烧反应延迟时间常数，T_F是燃油系统时间常数，T_CD是压缩排量时间常数。

步骤3.4：建立储能设备的低阶频率响应模型。储能设备的结构图如图8所示，其中，T_ESD为储能环节的储能时间常数。

步骤4：基于各资源低阶频率响应模型的分析后，为了提升计算速度，具体为整个系统扰动作为输入量，经聚合发电机惯量，输出为频差频差经各储能调节环节反馈至输入端。本发明建立了多资源、多机系统频率响应聚合模型。建模的核心思想是：先聚合多机，再聚合多资源。

步骤4.1：建立多机聚合的系统频率响应模型，火电机组的多机模型如图9所示，具体为整个系统扰动作为输入量，经聚合发电机惯量，输出为频差，频差经各机组调速器反馈至输入端。在多机系统频率响应模型中，每台机组的调差系数R_ti和机械功率增益K_ti与发电机的视在功率S_i密切相关，其中K_ti表示机组在整个电力系统发电机组的功率占比。因此，为了聚合模型，做了如下定义为：

其中,κ_ti为火电机组的等效增益参数。为了简化聚合模型，对调速器机组做归一化处理，定义每个火电机组支路的标准化增益λ_ti为：

F_H、T_R、T_G、T_C、T_R为聚合后频率响应模型的综合等效参数。在聚合计算中，需要找到火电参数最优解X＝{T_G,T_C,F_H,T_R},即聚合的频率响应模型参数可以等效为每一个机组参数的加权平均值。

水电机组的多机频率响应模型如图6所示，具体为整个系统扰动作为输入量，经聚合发电机惯量，输出为频差，频差经各机组调速器反馈至输入端。定义聚合模型的调差系数为：

其中,κ_hi为水电机组的等效增益参数。定义每个水电机组支路的标准化增益λ_hi为：

令燃气轮机参数最优解Y＝{T_w,T_y,T_RS,T_RH}，其系统参数聚合为：

燃气轮机组的多机频率响应模型如图7所示，具体为整个系统扰动作为输入量，经聚合发电机惯量，输出为频差，频差经各机组调速器反馈至输入端。定义聚合模型调差系数为：

其中,κ_gi为水电机组的等效增益参数。定义每个水电机组支路的标准化增益λ_gi为：

令水电参数最优解Z＝{c_g,b_g,X_g,Y_g,T_CR,T_F,T_CD}，其系统参数聚合为：

步骤4.2：将上述各资源模型多机分别聚合为单机模型后，进行多资源模型聚合。根据各资源模型的聚合理论，机械功率增益K与系统的容量占比密切相关，即：

其中，S_i为各资源机组的系统容量，S_sys为整个多资源机组的系统容量。

在发电机及负荷部分，以荷侧功率为基准，计算每种类型机组的容量占比为：

其中，H_i为各机组的发电机组的惯性时间常数，H为聚合后系统的发电机组的惯性时间常数。

在含有多资源、多机组电力系统频率响应聚合模型中，由于风电机组不具有实际的转子惯性，本发明不考虑其调频特性。依据“先聚合多机，再聚合多资源”的核心思想，构建多资源、多机系统频率响应聚合模型如图12所示，具体为：整个系统扰动作为输入量，经聚合发电机惯量，输出为频差，频差经各资源机组的聚合调速器反馈至输入端。

步骤5：为了提高计算速率，现对多资源、多机系统频率响应聚合模型进行开环解析。将上述含多种资源的电力系统频率响应聚合模型，进行了开环简化，并将各类型机组的调速器系统均拟合为一阶惯性环节。此外，本发明还忽略负荷阻尼D对频率的影响，其输出量为每台机组的有功出力，模型如图13所示，具体为：以系统负荷扰动为输入量，经发电机惯性环节输出为频差，将开环后频差信号经各资源机组的聚合调速器反馈至输入端。

根据开环模型，每台机组输出功率的计算如式所示：

式中，K_i表示各机组拟合后的功率—频率静特性系数；T_i表示各机组响应时间常数，ΔP_Gi为每台机组的调速器增发功率，Δω为系统的频差。在实际工程应用中，可认为在秒级尺度内系统频率初始衰减的斜率为常数m_Δω如图10所示，即：在扰动时间后，系统频差与时间成线性关系，解析式如下。

由此，对ΔP_G进行拉普拉斯逆变换可得其时域表达式分别为：

根据简化后的开环模型和叠加定理，Δω的时域表达式为：

令调速器拟合参数

积分后为：

当系统频率达到最低点时，有

此时时间为t_nadir，即频率最低点到达时间：

但上述式子中出现的C_i、t_nadir均为未知量，即共有N+1个未知量，故须建立如下的N+1次方程组进行求解。

最大频差Δω_max为：

频率最低点f_min为：f_min＝f₀-f_B·Δω_max (39)

式中，f₀为系统扰动前稳态频率，f_B为系统基础频率。

步骤6：将上述模拟计算出的系统失负荷量带入到系统频率响应模型中，求解系统的频率动态，其主要动态指标为：频率最大偏差、频率最低点到达时间以及频率的初始变化率。

步骤7：通过MATLAB编程方式实现上述过程，其仿真算例及其结果表明：本发明提出随机生产模拟方法考虑频率动态变化过程，与传统方法相比能够更加精确，能够满足面向快速频率响应的备用容量规划需求，为以低频减载为指标的备用容量规划提供了新思路。

具体实例，在火电机组的聚合模型验证中，以6台10MW机组为基准进行参数聚合，其结果见图15，取频率最低点，f_MM-SFR＝-0.3517Hz，f_ASFR＝-0.3555Hz。由此得出，Δf＝-0.0038Hz，其误差为：1.08％，该误差在可接受的范围内。在水电机组的聚合模型验证中，以2台50MW机组为基准进行参数聚合，其结果见图16，取频率最低点，f_MM-ASFR＝-0.4373Hz，f_ASFR＝-0.4321Hz。由此得出：Δf＝0.0052Hz，误差为：1.19％，在可接受的误差范围内。在燃气机组的聚合模型验证中，以1台150MW、1台100MW以及1台120MW机组为基准进行参数聚合，其结果见图17，取频率最低点，f_MM-ASFR＝-0.0246Hz，f_ASFR＝-0.0237Hz。由此得出：Δf＝0.0009Hz，误差为：3.7％，其误差在可接受的误差范围内。三组聚合模型算例以燃气轮机组为例，拟合调速器的一阶惯性环节如图18，其拟合误差在误差裕度内。在算例中，机组容量见表1，机组参数见表2。结合上述算例数据，模拟系统一年的失负荷量如图19所示，将准稳态频率作为参考指标与解析法和仿真法解出的频率最大频差对比，如图20所示。说明该方法提出的必要性和有效性，以及其解析法精度高，速度快的特点。

以上所述实施例仅表达本发明的实施方式，但并不能因此而理解为对本发明专利的范围的限制，应当指出，对于本领域的技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变形和改进，这些均属于本发明的保护范围。

Claims (5)

1.一种面向频率响应容量规划的电力系统随机生产模拟方法，其特征在于，所述模拟方法包括以下步骤：

步骤1：建立电力系统状态概率模型，模拟系统运行状况；所述电力系统状态概率模型包括风力发电机组的状态概率模型和常规发电机组的状态概率模型，具体步骤如下：

步骤1.1：根据马尔可夫链蒙特卡洛原理，建立风力发电机组的状态概率模型，包括风电机组运行状态转移与降额状态相互转换和运行状态与停运状态的相互转换，定义λ₁为风机系统由运行状态转移为停运状态的转移变化率，λ₂为风机系统由运行状态转移为降额状态的转移变化率，μ₁为风机系统由停运状态转移为运行状态的修复率，μ₂为风机系统由降额状态转移为运行状态的修复率；则根据风机的三状态模型以及马尔可夫理论，建立状态转移矩阵T：

基于马尔可夫过程的逼近原理，定义P为状态概率矩阵，I为单位矩阵且P＝[P₁,P₂,P₃]，P₁为风机处于运行状态的概率，P₂为风机处于降额状态的概率，P₃为风机处于停运状态的概率，得出：

对上式机型转置运算，并增加全概率条件，能够进一步得出风机各状态的概率：

在已知风机各状态概率的情况下，利用蒙特卡洛方法对系统进行随机抽样，即可确定风力发电系统的运行状态；

系统处在运行状态尚未进入降额状态的持续时间为τ₁₂：

系统处在运行状态尚未进入停运状态的持续时间为τ₁₃：

其中，MTTF为失效前平均持续运行时间，γ₁₂和γ₁₃为区间[0,1]内的均与分布随机数；

系统由降额状态或停运状态转移为正常运行状态的修复时间分别为τ₂₁和τ₃₁：

其中，MTTR为系统的平均修复时间，γ₂₁和γ₃₁为区间[0,1]内的均与分布随机数；

步骤1.2：同样，根据马尔可夫链蒙特卡洛原理，建立常规机组的状态概率模型，主要为运行转态和停运状态之间的转换，定义FOR为强迫停运率，即当蒙特卡洛随机抽样结果小于该值，系统进入停运状态；当其结果大于该值时，系统处于正常运行状态；其求解方法如下：

其系统在每个状态下的运行时间和风力发电机组的算法相似，正常运行状态下和停运状态下的运行时间分别为τ₁和τ₂，即：

τ₁＝-MTTF lnγ₁ (11)

τ₂＝-MTTR lnγ₂ (12)

其中，γ₁和γ₂为区间[0,1]之间均匀分布的随机数；

步骤1.3:引入故障过滤的算法，解决蒙特卡洛法在进行随机抽样的过程中出现连续停运的问题；

步骤2：通过状态概率模型，得出各类型机组的生产状态，结合其工况建立电力系统输出功率模型，模拟系统的失负荷量；包括建立风力发电输出功率模型、建立常规机组输出功率模型以及建立失负荷量模型；具体为：

步骤2.1：基于ARMA模型建立风力发电输出功率模型；

基于风速的历史数据，对其进行标准化，得出风速平均值μ和方差σ，t时刻的预测风速SW_t的表达式如下：

SW_t＝μ+σy_t (13)

其中，y_t为预测的时间序列值；

假设系统扰动的白噪声序列{ε_t}的均值为0，方差为

则时间序列y_t的为：

其中，

为自回归系数，θ_i(j＝1,2,3,…,q)为滑动平均系数，p和q分别ARMA模型的自回归阶数和滑动平均阶数；

采用AIC准则对ARMA模型进行模型的阶数确定；

在AIC准则的模型中，当p和q达到某一数值时，使得AIC的值最小；此时，使得AIC的值达到最小的p和q为ARMA模型的阶数；

在风力发电机组中，定义v_ci是时风力机组启动时的切入风速，v_r是风力机组正常运行时的额定风速，v_co是风力机组运行时的切出风速，保证风机安全使风机停止工作；P_r是额定风速下风力发电机组正常运行时的额定功率；当风速超过切入风速v_ci且未到达额定风速v_r时，风机的输出功率与风速是非线性关系；当风速超过额定风速v_r但未到达切出风速v_co时，风机保持额定功率一直不变；当风速超过切出风速v_co时，风机停止工作；将上述过程建立为如下的分段函数：

其中，P_w为t时刻风电机组的输出功率，v_t为t时刻风机的预测风速，α、β、γ为风电机组输出功率的待定系数，

步骤2.2：建立常规机组的输出功率模型，其中常规发电机组包括火力发电机组、水力发电机组和燃气轮机发电机组，以算例IEEE RTS-79中各机组的系统容量为基准，计算其输出功率；

步骤2.3在电力系统的总输出功率的部分，包括常规和风电机组两个部分，将系统的输出功率与系统状态相结合，分别进行计算，即：

其中，P为系统的总输出功率，P_i为各类型机组的输出功率，α_i为各类型机组的系统状态；

因此，电力系统的失负荷量ΔP_d为：

ΔP_d＝P-P_L (19)

其中，P_L为算例中的电力系统负荷；

步骤3：将上述模拟计算出的系统失负荷量带入系统频率响应模型中，求解系统的频率动态；因此，需建立各资源下系统的频率响应模型；从低阶频率响应模型出发，进行模型的建立与分析；

步骤3.1：建立火电机组的低阶频率响应模型；

机组功率增量高压缸汽轮机惯性环节和机械功率增益环节，与负荷侧扰动一同输入到发电机惯性环节中，输出频差，最后频差经调速器负反馈环节，反馈到输入端；其中，定义P_e为发电机的负荷侧功率、P_m为汽轮机的机械功率、P_SP为机组功率增量，P_α为系统加速功率；T_R为再热时间常数、H为系统惯性常数、1/R_T负反馈增益、R_T为调差系数、F_H为高压汽轮机功率比例系数、K_T为机械功率增益、D为阻尼系数、T₁为调速器时间常数、T₂为汽容时间常数；

步骤3.2：建立水电机组的低阶频率响应模型；

负荷侧扰动，经水电机组的发电机转子惯性环节，输出频差，最后频差经调速器负反馈环节至输入端；定义H为机组的惯性时间常数，D为机组的阻尼常数，T_y为接力器的响应时间常数，T_w为水流的惯性时间常数，T_RH为水轮机调速器暂态下垂时间常数，T_RS为水轮机调速器复位时间，1/R_h为水轮机的调差系数；

步骤3.3：建立燃气轮机组的低阶频率响应模型；

负荷侧扰动，经燃气轮机机组的发电机转子惯性环节，输出频差，最后频差经调速器负反馈环节至输入端；其中，定义R_G为系统的调差系数，K_G为功率容量占比，c_g和b_g为燃油系统中阀门定位器的传递函数的系数，X_g是调速器部分的超前时间常数，Y_g是调速器部分的滞后时间常数，T_CR是燃烧反应延迟时间常数，T_F是燃油系统时间常数，T_CD是压缩排量时间常数；

步骤3.4：建立储能设备的低阶频率响应模型，T_ESD为储能环节的储能时间常数；

步骤4：基于各资源低阶频率响应模型的分析后，为提升计算速度，为整个系统扰动作为输入量，经聚合发电机惯量，输出为频差频差经各储能调节环节反馈至输入端；通过先聚合多机，再聚合多资源的方式建立多资源、多机系统频率响应聚合模型；具体如下：

步骤4.1：建立多机聚合的系统频率响应模型

整个系统扰动作为输入量，经聚合发电机惯量，输出为频差，频差经各机组调速器反馈至输入端；在多机系统频率响应模型中，每台机组的调差系数R_ti和机械功率增益K_ti与发电机的视在功率S_i密切相关，其中K_ti表示机组在整个电力系统发电机组的功率占比；因此，为了聚合模型，做了如下定义为：

其中,κ_ti为火电机组的等效增益参数；为简化聚合模型，对调速器机组做归一化处理，定义每个火电机组支路的标准化增益λ_ti为：

其中,F_H、T_R、T_G、T_C、T_R为聚合后频率响应模型的综合等效参数；聚合的频率响应模型参数可以等效为每一个机组参数的加权平均值；

水电机组的多机频率响应模型为整个系统扰动作为输入量，经聚合发电机惯量，输出为频差，频差经各机组调速器反馈至输入端；定义聚合模型的调差系数为：

其中,κ_hi为水电机组的等效增益参数；定义每个水电机组支路的标准化增益λ_hi为：

令燃气轮机参数最优解Y＝{T_w,T_y,T_RS,T_RH}，其系统参数聚合为：

燃气轮机组的多机频率响应模型为整个系统扰动作为输入量，经聚合发电机惯量，输出为频差，频差经各机组调速器反馈至输入端；定义聚合模型调差系数为：

其中,κ_gi为水电机组的等效增益参数；定义每个水电机组支路的标准化增益λ_gi为：

令水电参数最优解Z＝{c_g,b_g,X_g,Y_g,T_CR,T_F,T_CD}，其系统参数聚合为：

步骤4.2：将上述各资源模型多机分别聚合为单机模型后，进行多资源模型聚合；根据各资源模型的聚合理论，得到：

其中，S_i为各资源机组的系统容量，S_sys为整个多资源机组的系统容量；

在发电机及负荷部分，以荷侧功率为基准，计算每种类型机组的容量占比为：

其中，H_i为各机组的发电机组的惯性时间常数，H为聚合后系统的发电机组的惯性时间常数；

在含有多资源、多机组电力系统频率响应聚合模型中，不考虑其调频特性；依据“先聚合多机，再聚合多资源”的方式构建多资源、多机系统频率响应聚合模型，具体为：整个系统扰动作为输入量，经聚合发电机惯量，输出为频差，频差经各资源机组的聚合调速器反馈至输入端；

步骤5：对多资源、多机系统频率响应聚合模型进行开环解析；将上述含多种资源的电力系统频率响应聚合模型，进行开环简化，并将各类型机组的调速器系统均拟合为一阶惯性环节；此外，忽略负荷阻尼D对频率的影响，其输出量为每台机组的有功出力，模型为：以系统负荷扰动为输入量，经发电机惯性环节输出为频差，将开环后频差信号经各资源机组的聚合调速器反馈至输入端；

根据开环模型，每台机组输出功率的计算如式所示；

式中，K_i表示各机组拟合后的功率—频率静特性系数；T_i表示各机组响应时间常数，ΔP_Gi为每台机组的调速器增发功率，Δω为系统的频差；

当系统频率达到最低点时，此时时间为t_nadir，即频率最低点到达时间：

建立如下的N+1次方程组对上式进行求解；

最大频差Δω_max为：

频率最低点f_min为：

f_min＝f₀-f_B·Δω_max (39)

式中，f₀为系统扰动前稳态频率，f_B为系统基础频率；

步骤6：将上述模拟计算出的系统失负荷量带入到系统频率响应模型中，求解系统的频率动态，其主要动态指标为：频率最大偏差、频率最低点到达时间以及频率的初始变化率；

步骤7：通过MATLAB方式实现上述过程。

2.根据权利要求1所述的一种面向频率响应容量规划的电力系统随机生产模拟方法，其特征在于，含有多资源、多机组的系统频率响应聚合模型，在模型中以单机的形式存在，主要包括：风电、火电、水电、燃气轮机以及储能设备。

3.根据权利要求1所述的一种面向频率响应容量规划的电力系统随机生产模拟方法，其特征在于，步骤1所述步骤1.3中，引入故障过滤算法具体过程为：

1)第i台机组进入故障过滤检测，检测其在整个周期内，第j个时刻是否有连续故障趋势；若没有，则跳出故障过滤程序；

2)若具有连续故障趋势，搜寻从该时刻开始共有多少个故障时刻；若超出规定时间，结束搜索计算，将该时刻之后，规定时刻之前的状态赋值0；若未超出规定时间，将该时刻之后的所有故障连续时刻状态赋值0；

3)继续检测故障时刻，直到每台机组的所有时刻的故障状态全部检测完毕，结束故障过滤。

4.根据权利要求1所述的一种面向频率响应容量规划的电力系统随机生产模拟方法，其特征在于，步骤2所述步骤2.1中AIC准则的函数表达如下：

AIC(p,q)＝logσ²(p,q)+2(p+q+1)/N (15)

式中，σ²(p,q)为预测时间序列的残差方差，N为样本的个数。

5.根据权利要求1所述的一种面向频率响应容量规划的电力系统随机生产模拟方法，其特征在于，步骤2.1中公式(19)中待定系数的求解算式如下：

Images