CN107846024B - 孤岛微电网概率潮流解析计算方法 - Google Patents

Abstract

本发明提供一种孤岛微电网概率潮流的解析计算方法。为包含光伏、风电等间歇性DG以及负荷随机波动的孤岛微电网的运行状态分析提供计算方法。现有的孤岛微电网稳态潮流分析方法主要集中在常规潮流计算，仅能获取单一运行点的潮流结果，而关于概率潮流计算的研究也仅限于模拟法，即需要根据大量抽样样本计算统计获得概率潮流结果，计算时间长不利于在线分析，且多数未考虑微电网中DG出力的相关性。本发明方法针对这些局限，提供一种孤岛微电网概率潮流的解析计算方法，能够实现孤岛微电网概率潮流的快速计算，且精度能够符合工程精度要求，具备在线分析应用的可行性。

Description

孤岛微电网概率潮流解析计算方法

技术领域

本发明属于孤岛微电网稳态分析领域，尤其涉及一种孤岛微电网概率潮流解析计算方法。

背景技术

微电网(microgrid，MG)为综合利用可再生能源提供了技术手段，根据与主网有无电气联系又可分为并网型与孤岛型，其中孤岛型微电网由于无平衡节点，缺乏主网支撑，通常采用对等控制方式，即由各采取下垂控制策略的分布式电源 (distributed generator，DG)、储能装置(energy storage，ES)、静止无功发生器(STATCOM)等共同参与微电网的频率与电压的调节。

孤岛微电网潮流计算是对孤岛微电网规划、运行分析、优化调度的基础，目前已有专利或文献建立了孤岛微电网的潮流模型，并通过信赖域算法、类奔德斯分解方法等提高潮流计算收敛性与收敛速度。诸如光伏、风电这一类DG的出力具有较强的随机性，同时负荷也存在波动，相比于单一运行点的潮流结果，概率潮流(probabilistic load flow，PLF)能够获得各状态变量的概率分布信息，更为全面地反映微电网的运行状态。

求解PLF的方法主要分为模拟法及解析法：模拟法主要是指蒙特卡洛法，其缺陷是计算量过大，通常用于校验概率潮流准确性，不宜在线计算或用于优化调度。而解析法则是基于近似线性化的数学假设，根据不确定量数值特征求取近似解析表达式，仅需要在基准点进行一次潮流计算，计算效率高。目前有文献基于半不变量法求解主从控制模式下的孤岛微电网PLF，计算时将柴油发电机视为平衡节点，储能视为恒功率负荷，与传统电力系统概率潮流计算类似，但未考虑对等控制方式下孤岛微电网的运行特点。

发明内容

针对现有技术存在的不足和空白，本发明采用以下技术方案：

一种孤岛微电网概率潮流解析计算方法，其特征在于，包括以下步骤：

步骤一：输入原始数据，包括：微电网网架数据、电源和负荷的基本参数以及所对应的概率模型；

步骤二：由已知m个输入随机变量W＝[w₁,w₂,……,w_m]的边际概率分布F(W)及相关系数矩阵R_W建立满足给定相关性水平的输入变量概率模型；

步骤三：根据各节点DG出力概率模型与负荷概率模型计算各节点功率注入期望值，基于孤岛微电网潮流模型进行基准运行点的稳态潮流计算，获得基准点潮流结果；

步骤四：计算各输入随机变量的各阶半不变量，并将具有相关性的随机变量的各阶半不变量转换为独立随机变量的各阶半不变量；

步骤五：计算节点状态变量、支路潮流变量以及具有下垂控制特性的电源输出功率的各阶半不变量；

步骤六：由Gram-Charlier级数展开获取输出变量的概率密度函数(probabilitydensity function，PDF)，及累积分布函数(cumulative distribution function， CDF)，输出PLF计算结果。

进一步地，步骤二的具体步骤为：

步骤二A：已知m个输入随机变量W＝[w₁,w₂,……,w_m]的边际概率分布F(W)，并使用相关系数矩阵R_W描述出力相关性，所述相关系数矩阵表达式为：

矩阵中各元素按下式计算：

(2)式中：

为随机变量w_i与w_j之间的相关系数；Cov(w_i,w_j)为w_i与w_j的协方差；σ_wi与σ_wj为w_i与w_j的标准差；

建立满足给定相关性水平的输入变量概率模型的步骤为：

步骤二B：产生满足独立正态分布的样本矩阵S，其维数为m×N，N为样本规模；步骤二C：运用三阶多项式正态变换(third-order polynomial normal transformation，TPNT)以及求解代数方程的方法求取标准正态分布随机变量Y 的等效相关系数矩阵R_Y；

步骤二D：对R_Y进行Cholesky分解获得其下三角矩阵L_Y，并由Y＝L_YS获得具有相关性的服从标准正态分布的样本矩阵Y；

步骤二E：通过等概率转换原则得到满足相关性水平的样本矩阵W，从而获得满足给定相关性水平的输入变量概率模型。

进一步地，步骤三中，所述基于孤岛微电网潮流模型进行基准运行点的稳态潮流计算包括以下步骤：

步骤三A：设定节点功率方程统一表达式：仅设定PQ节点与PV节点两种节点类型，其中PQ节点需列写有功与无功平衡方程，PV节点仅需列写有功平衡方程，方程形式为：

式中：P_Ci和Q_Ci分别为节点i恒功率电源注入的有功和无功；P_Di和Q_Di分别为节点i具有下垂控制特性的设备注入的有功和无功；P_Li和Q_Li分别为节点i的有功与无功负荷；P_i和Q_i分别为节点i注入的总有功功率和无功功率；

其中，所述具有下垂控制特性的设备注入的功率可统一表示为：

(4)、(5)、(6)、(7)式中：P_Dimax、P_Dimin、Q_Dimax、Q_Dimin分别为节点i具有下垂控制特性的设备注入的有功的上限、节点i具有下垂控制特性的设备注入的有功的下限、节点i具有下垂控制特性的设备注入的无功的上限和节点i具有下垂控制特性的设备注入的无功的下限，若无具有下垂控制特性的设备则均取为0； f_max、f_min、U_max、U_min分别为系统频率允许的上限、系统频率允许的下限、电压允许的上限、电压允许的下限；K_Dfi、K_DUi分别为节点i对应的P-f、Q-U下垂系数；f₀、f、U₀、U_i分别为系统频率的空载值、系统频率、电压的空载值与电压在节点i的实际值；

其中，对负荷进行建模可表示为：

(8)式中：P_LNi和Q_LNi分别为节点i在额定工况下的有功与无功负荷；U_Ni和f_N分别为节点i的额定电压与频率，取U_Ni＝1p.u.，f_N＝1p.u.；P_i和Q_i分别为节点i注入的总有功功率和无功功率；A_pi、B_pi、C_pi、A_qi、B_qi、C_qi分别为节点i的负荷有功功率中恒阻抗型的百分比系数、节点i的负荷有功功率中恒电流型的百分比系数、节点i的负荷有功功率中恒功率型的百分比系数、节点i的负荷无功功率中恒阻抗型的百分比系数、节点i的负荷无功功率中恒电流型的百分比系数、节点i的负荷无功功率中恒功率型的百分比系数，各满足总和为1；k_Lpi、 k_Lqi分别为节点i负荷的有功和无功功率的静态频率调节系数；U_i为节点i系统电压的实际值；f为节点i系统频率的实际值；

其中，节点注入的有功与无功功率可表示为：

(9)式中：n为节点数目；G_ij与B_ij分别为节点导纳矩阵的实部与虚部；δ_ij为节点i和节点j的相角差；U_i为节点i系统电压的实际值；U_j为节点j系统电压的实际值；

步骤三B：牛顿拉夫逊法求解:孤岛微电网中PQ节点与PV节点的数目分别为n_PQ与n_PV，则总计列写2n_PQ+n_PV个方程，采用牛顿拉夫逊法求解潮流方程组，其中，修正方程简写为：

(10)式中：ΔP、ΔQ为节点有功与无功不平衡量；Δf、Δδ、ΔU为分别为频率、相角、电压的修正量；J为雅克比矩阵，其分块矩阵分别为：

雅克比矩阵J中各分块矩阵的维数为：

各分块矩阵的表达式为：

(13)式中，各参数由(3)-(9)式中各参数代入计算。

进一步地，步骤四、步骤五的计算步骤为：

步骤四五A：在基准运行点处将潮流方程线性化并忽略2次及以上高次项，可得随机分量的关系为：

(14)式中：ΔW、ΔX、ΔZ分别为节点注入功率、节点状态变量以及支路潮流变量的随机分量；J₀为基准点处的雅克比矩阵；X₀为基准点处的节点状态量；

当已知W的各阶半不变量，则可根据半不变量的线性性质获得X、Z的各阶半不变量；

步骤四五B：通过对R_W进行Cholesky分解获得其下三角矩阵L_W，可将m个相关输入随机变量W转换为m个相互独立的变量W_ind的组合，转换关系如下：

W_ind＝L_W ^-1W (15)

相应的J₀中与原始m个相关输入随机变量W对应列的元素需要通过L_W修正；

步骤四五C：具有下垂控制特性的电源输出功率与系统频率及节点电压有关，为获取其输出的概率分布，需要计算具有下垂控制特性的设备输出功率的半不变量，假定在系统运行范围内，对式(4)(6)在基准点线性化，可推得其关系为：

结合式(14)进一步转换为与节点注入功率ΔW的关系后同样可基于半不变量的线性性质由W的各阶半不变量求得具有下垂控制特性的电源输出功率的各阶半不变量。

进一步地，步骤六中获取输出变量的概率密度函数的步骤为：

步骤六A：借助标准正态分布概率密度函数φ以及它的各阶导数的线性组合来逼近最终待求的分布函数，根据已求得状态变量的各阶半不变量为γ_k，列出由各阶半不变量推导各阶中心矩的公式：

步骤六B：由半不变量γ_k求得中心矩后，进一步通过中心矩表示Gram-Charlier级数系数，公式为：

步骤六C：将随机变量标准化处理后借助Gram-Charlier级数展开可将随机变量的分布函数表达为由正态随机变量各阶导数组成的级数；

则由Gram-Charlier级数展开获取输出变量的概率密度函数及累积分布函数，输出PLF计算结果。

本发明的主要贡献是提出一种对等控制模式下孤岛微电网概率潮流的解析计算方法。首先综合DG、ES、STATCOM等设备的下垂控制特性建立节点功率方程统一表达式，简化运算与编程；随后在孤岛微电网潮流模型的基础上，扩展出概率潮流模型，运用Nataf变换对微电网中DG的出力相关性进行处理；最后基于半不变量与Gram-Charlier级数展开获得系统频率及节点电压等状态变量的完整概率分布信息。

本发明的主要目的提供一种孤岛微电网概率潮流的解析计算方法。为包含光伏、风电等间歇性DG以及负荷随机波动的孤岛微电网的运行状态分析提供计算方法。现有的孤岛微电网稳态潮流分析方法主要集中在常规潮流计算，仅能获取单一运行点的潮流结果，而关于概率潮流计算的研究也仅限于模拟法，即需要根据大量抽样样本计算统计获得概率潮流结果，计算时间长不利于在线分析，且多数未考虑微电网中DG出力的相关性。本发明方法针对这些局限，提供一种孤岛微电网概率潮流的解析计算方法，能够实现孤岛微电网概率潮流的快速计算，且精度能够符合工程精度要求，具备在线分析应用的可行性。

本发明有益效果包括：

1)综合考虑DG、ES、STATCOM等设备的下垂控制特性建立节点功率方程统一表达式，简化了节点分类与编程，在系统结构不发生变化情况下，即使发生设备投退的情况，也可不必修改潮流方程与雅克比矩阵结构，只需更新数值。

2)本发明方法的计算结果与蒙特卡洛法计算结果接近，同时在DG容量、出力相关性以及负荷波动的变化下，误差均在可接受范围内，能够满足工程精度要求。

3)本发明方法避免了抽样计算，计算时间短、效率高，与模拟法相比具有显著优势，具备应用于在线计算的可行性。因此在孤岛微电网运行分析、风险预警、优化调度等领域具有应用的可行性。

附图说明

下面结合附图和具体实施方式对本发明进一步详细的说明：

图1为本发明实施例流程示意图；

图2为本发明实施例中算例系统网络节点示意图；

图3为本发明实施例中孤岛微电网频率PDF与CDF统计图；

图4为本发明实施例中节点1～18的电压期望值与置信区间统计图；

图5为本发明实施例中节点9STATCOM输出无功的PDF及CDF统计图；

图6为本发明实施例中节点18储能装置输出有功的PDF及CDF统计图；

图7为本发明实施例中调整前后的频率PDF与柴油发电机负载率PDF对比图。

具体实施方式

为让本专利的特征和优点能更明显易懂，下文特举实施例，作详细说明如下：

如图1所示，本发明实施例包括以下步骤：

步骤一：输入原始数据，包括：微电网网架数据、电源和负荷的基本参数以及所对应的概率模型；

步骤二：由已知m个输入随机变量W＝[w₁,w₂,……,w_m]的边际概率分布F(W)及相关系数矩阵R_W建立满足给定相关性水平的输入变量概率模型；

步骤三：根据各节点DG出力概率模型与负荷概率模型计算各节点功率注入期望值，基于孤岛微电网潮流模型进行基准运行点的稳态潮流计算，获得基准点潮流结果；

步骤四：计算各输入随机变量的各阶半不变量，并将具有相关性的随机变量的各阶半不变量转换为独立随机变量的各阶半不变量；

步骤五：计算节点状态变量、支路潮流变量以及具有下垂控制特性的电源输出功率的各阶半不变量；

步骤六：由Gram-Charlier级数展开获取输出变量的概率密度函数(probabilitydensity function，PDF)，及累积分布函数(cumulative distribution function， CDF)，输出PLF计算结果。

其中，步骤二的具体步骤为：

步骤二A：已知m个输入随机变量W＝[w₁,w₂,……,w_m]的边际概率分布F(W)，并使用相关系数矩阵R_W描述出力相关性，相关系数矩阵表达式为：

矩阵中各元素按下式计算：

(2)式中：

为随机变量w_i与w_j之间的相关系数；Cov(w_i,w_j)为w_i与w_j的协方差；σ_wi与σ_wj为w_i与w_j的标准差；

建立满足给定相关性水平的输入变量概率模型的步骤为：

步骤二B：产生满足独立正态分布的样本矩阵S，其维数为m×N，N为样本规模；步骤二C：运用三阶多项式正态变换(third-order polynomial normal transformation，TPNT)以及求解代数方程的方法求取标准正态分布随机变量Y 的等效相关系数矩阵R_Y；

步骤二D：对R_Y进行Cholesky分解获得其下三角矩阵L_Y，并由Y＝L_YS获得具有相关性的服从标准正态分布的样本矩阵Y；

步骤二E：通过等概率转换原则得到满足相关性水平的样本矩阵W，从而获得满足给定相关性水平的输入变量概率模型。

其中，步骤三中，基于孤岛微电网潮流模型进行基准运行点的稳态潮流计算包括以下步骤：

步骤三A：设定节点功率方程统一表达式：仅设定PQ节点与PV节点两种节点类型，其中PQ节点需列写有功与无功平衡方程，PV节点仅需列写有功平衡方程，方程形式为：

式中：P_Ci和Q_Ci分别为节点i恒功率电源注入的有功和无功，例如光伏、风电等采取MPPT控制方式的DG属于此类电源；P_Di和Q_Di分别为节点i具有下垂控制特性的设备注入的有功和无功，例如柴油发电机、燃气轮机、储能、STATCOM等均属于此类电源；P_Li和Q_Li分别为节点i的有功与无功负荷；P_i和Q_i分别为节点i 注入的总有功功率和无功功率；

其中，具有下垂控制特性的设备注入的功率可统一表示为：

(4)、(5)、(6)、(7)式中：P_Dimax、P_Dimin、Q_Dimax、Q_Dimin分别为节点i具有下垂控制特性的设备注入的有功和无功的上下限，不同设备所对应的参数取值及含义见表1，若无具有下垂控制特性的设备则均取为0；f_max、f_min、U_max、U_min分别为系统频率、电压允许上下限；K_Dfi、K_DUi分别为节点i对应的P-f、Q-U 下垂系数；f₀、f、U₀、U_i分别为系统频率、电压的空载值与节点i的实际值；

其中储能可以工作在充放电两种工况，在不考虑荷电状态的情况下，其 P_Dimin、P_Dimax取决于储能的最大充放电功率，同样STATCOM可以吸收或发出无功， Q_Dimin、Q_Dimax取决于STATCOM配置容量，根据方程符号及物理含义，充电功率最大值以及STATCOM吸收无功最大值为一负值。

其中，对负荷进行建模时计及负荷的电压和频率静特性，负荷使用恒阻抗、恒电流和恒功率的组合模型来描述，可表示为：

(8)式中：P_LNi和Q_LNi分别为节点i在额定工况下的有功与无功负荷；U_Ni和f_N分别为节点i的额定电压与频率，取U_Ni＝1p.u.，f_N＝1p.u.；P_i和Q_i分别为节点i注入的总有功功率和无功功率；A_pi、B_pi、C_pi、A_qi、B_qi、C_qi分别为节点 i的负荷有功与无功功率中恒阻抗型、恒电流型、恒功率型的百分比系数，各满足总和为1；k_Lpi、k_Lqi分别为节点i负荷的有功和无功功率的静态频率调节系数； U_i为节点i系统电压的实际值；f为节点i系统频率的实际值；

其中，节点注入的有功与无功功率可表示为：

(9)式中：n为节点数目；G_ij与B_ij分别为节点导纳矩阵的实部与虚部；δ_ij为节点i和节点j的相角差；U_i为节点i系统电压的实际值；U_j为节点j系统电压的实际值；

步骤三B：牛顿拉夫逊法求解:孤岛微电网中PQ节点与PV节点的数目分别为n_PQ与n_PV，则总计列写2n_PQ+n_PV个方程，采用牛顿拉夫逊法求解潮流方程组，其中，修正方程简写为：

(10)式中：ΔP、ΔQ为节点有功与无功不平衡量；Δf、Δδ、ΔU为分别为频率、相角、电压的修正量；J为雅克比矩阵，其分块矩阵分别为：

雅克比矩阵J中各分块矩阵的维数为：

各分块矩阵的表达式为：

(13)式中，各参数由(3)-(9)式中各参数代入计算。

其中，步骤四、步骤五的计算步骤为：

步骤四五A：在基准运行点处将潮流方程线性化并忽略2次及以上高次项，可得随机分量的关系为：

(14)式中：ΔW、ΔX、ΔZ分别为节点注入功率、节点状态变量以及支路潮流变量的随机分量；J₀为基准点处的雅克比矩阵；X₀为基准点处的节点状态量；

当已知W的各阶半不变量，则可根据半不变量的线性性质获得X、Z的各阶半不变量；

步骤四五B：通过对R_W进行Cholesky分解获得其下三角矩阵L_W，可将m个相关输入随机变量W转换为m个相互独立的变量W_ind的组合，转换关系如下：

W_ind＝L_W ^-1W (15)

相应的J₀中与原始m个相关输入随机变量W对应列的元素需要通过L_W修正；

步骤四五C：具有下垂控制特性的电源输出功率与系统频率及节点电压有关，为获取其输出的概率分布，需要计算具有下垂控制特性的设备输出功率的半不变量，假定在系统运行范围内，对式(4)(6)在基准点线性化，可推得其关系为：

结合式(14)进一步转换为与节点注入功率ΔW的关系后同样可基于半不变量的线性性质由W的各阶半不变量求得具有下垂控制特性的电源输出功率的各阶半不变量。

其中，步骤六中获取输出变量的概率密度函数的步骤为：

步骤六A：借助标准正态分布概率密度函数φ以及它的各阶导数的线性组合来逼近最终待求的分布函数，根据已求得状态变量的各阶半不变量为γ_k，由于 Gram-Charlier级数的各项系数与中心矩满足一定的关系，因此先列出由各阶半不变量推导各阶中心矩的公式：

步骤六B：由半不变量γ_k求得中心矩后，进一步通过中心矩表示Gram-Charlier级数系数，公式为：

步骤六C：将随机变量标准化处理后借助Gram-Charlier级数展开可将随机变量的分布函数表达为由正态随机变量各阶导数组成的级数。

则由Gram-Charlier级数展开获取输出变量的概率密度函数及累积分布函数，输出PLF计算结果。

具体地，本发明实施例运行在配置为IntelPentium双核处理器，2.0G内存的计算机上，于MATLABR2013a环境下编写程序实现本发明所提出的方法。以 IEEE33节点系统作为算例，将节点1与主网断开并配置DG，形成孤岛微电网。

算例系统网络节点如图2所示。将系统基准容量取1MVA，取孤岛微电网的安全运行范围为f_max＝1.004p.u.，f_min＝0.996p.u.，U_max＝1.05p.u.，U_min＝0.95p.u.。各节点恒阻抗、恒电流、恒功率负荷占比统一取为0.3、0.3、0.4，静态频率调节系数取k_Lpi＝2、k_Lqi＝-2。

负荷服从正态分布，负荷期望值见表1，负荷波动的标准差σ_load取为期望值的10％。

表1算例系统负荷数据

在节点1、9、18接入同等容量的光伏电源，假定待分析时段内各光伏电源出力均满足在区间[400kW，600kW]范围的Beta分布，形状参数为α＝1.693，β＝5.162，各光伏电源之间的相关系数ρ_DG均取为0.5。其他接入设备的容量参数见表2，各设备均采取下垂控制策略，下垂系数根据式(3)(5)进行整定。

表2设备参数

以蒙特卡洛法进行10000次抽样计算所得结果作为准确值，分析本发明实施例方法的误差情况，各节点电压与系统频率的误差情况见表3。

表3各节点电压与系统频率的误差情况

结果表明本发明实施例方法与蒙特卡洛法结果接近，能够满足工程精度的要求。利用本发明实施例方法以及蒙特卡洛法统计求出的孤岛微电网频率PDF与 CDF如图3所示。

本发明实施例方法计算得PDF与CDF和蒙特卡洛法统计得到的结果十分接近，系统基本运行在额定频率附近，受光伏及负荷波动性影响可能存在偏离额定频率的情况，但偏离程度不大。

通过PLF计算可获得各节点电压的概率分布信息，节点1～18的电压期望值与置信区间统计图参见图3。

由图4可知，在光伏及负荷存在波动的环境下，PLF能够直观显示孤岛微电网可能得运行范围，并给出具有概率意义的信息，可为孤岛微电网的运行提供参考与决策依据。

此外，验证本发明实施例方法计算采取下垂控制策略的设备输出情况的准确性，其中节点9的STATCOM输出无功以及节点18储能装置输出有功的PDF与CDF 参见图5与图6。

由图5和图6可知，本发明实施例方法所获得的设备输出功率的PDF与CDF 与蒙特卡洛法所得结果基本接近。考虑节点9的STATCOM配置容量，由于并网点电压低于额定值，因此工作在发出无功的状态，且伴随着节点9电压的不确定性，同样输出无功也呈现不确定性。

而节点18的储能装置基本工作在以较小功率充电或放电的状态，主要由于此时系统频率基本在额定值附近波动，且波动范围不大，相对而言储能也在充放电之间波动，且功率同样较小。由此可见本发明实施例方法能够较为准确地估计出具有下垂控制特性设备的可能输出范围，为评估设备运行状态提供参考。

在原算例条件的基础上，分别单独改变算例中光伏装机容量S_DG、出力相关性ρ_DG以及负荷波动标准差σ_load，分析不同因素对算法误差的影响。其中表中 S_DG的数值指相对于原算例的光伏装机比例，表中ε_μmean、ε_μmax表示期望值的平均与最大偏差绝对值，ε_σmean、ε_σmax表示标准差的平均与最大偏差绝对值。不同场景下的误差变化情况见表4。

表4不同场景下的误差情况

可见随着并网光伏容量增加、出力相关性增强以及负荷波动增强，相应的算法误差都会有所增大，主要是由于半不变量法在基准运行点线性化，当输入变量的不确定性增强会引起误差增大。但数据表明误差均处于可接受范围内，都能够满足工程精度要求。

而在计算时间上，本发明实施例方法仅需0.06s即可完成PLF计算，而抽样规模为10000次的蒙特卡洛模拟法需要163.27s，因此本发明实施例方法在允许一定误差的前提下能够在计算时间上具有显著优势，具备应用于在线计算的可行性。

进一步针对上述算例系统就本发明PLF方法在孤岛微电网的运行分析中的应用进行举例说明。以柴油发电机组的启停为例，原算例中节点33有4台额定功率500kW的柴油发电机并网运行，系统频率与柴油发电机负载率的PDF如图中浅色曲线所示。柴油发电机最低负载要求通常在30％左右，当负载率较低时发电成本较高，经济效益差。由图7可知，柴油发电机此时负载率大致处于30％～45％的范围，负载率较低，考虑进行调整，运用本发明实施例方法计算退出1台柴油发电机后的PLF，调整后系统频率与柴油发电机负载率PDF见深色曲线。

通过对比可见，尽管系统频率略低于调整前，但仍大致处于0.999 p.u.～1.000p.u.之间，频率完全合格，而柴油发电机的负载率均提升到55％左右，较为接近最佳负载状态，经济性较好。通过PLF结果可以判定调整后在系统安全性不受到影响的情况下能够提高经济性，因此可以令1台柴油发电机退出运行。

本专利不局限于上述最佳实施方式，任何人在本专利的启示下都可以得出其它各种形式的孤岛微电网概率潮流解析计算方法，凡依本发明申请专利范围所做的均等变化与修饰，皆应属本专利的涵盖范围。

Claims (5)

1.一种孤岛微电网概率潮流解析计算方法，其特征在于，包括以下步骤：

步骤一：输入原始数据，包括：微电网网架数据、电源和负荷的基本参数以及所对应的概率模型；

步骤二：由已知m个输入随机变量W＝[w₁,w₂,……,w_m]的边际概率分布F(W)及相关系数矩阵R_W建立满足给定相关性水平的输入变量概率模型；

步骤三：根据各节点DG出力概率模型与负荷概率模型计算各节点功率注入期望值，基于孤岛微电网潮流模型进行基准运行点的稳态潮流计算，获得基准点潮流结果；

步骤四：计算各输入随机变量的各阶半不变量，并将具有相关性的随机变量的各阶半不变量转换为独立随机变量的各阶半不变量；

步骤五：计算节点状态变量、支路潮流变量以及具有下垂控制特性的电源输出功率的各阶半不变量；

步骤六：由Gram-Charlier级数展开获取输出变量的概率密度函数(probabilitydensity function，PDF)，及累积分布函数(cumulative distribution function，CDF)，输出概率潮流计算结果。

2.根据权利要求1所述的孤岛微电网概率潮流解析计算方法，其特征在于，步骤二的具体步骤为：

步骤二A：已知m个输入随机变量W＝[w₁,w₂,……,w_m]的边际概率分布F(W)，并使用相关系数矩阵R_W描述出力相关性，所述相关系数矩阵表达式为：

矩阵中各元素按下式计算：

(2)式中：

为随机变量w_i与w_j之间的相关系数；Cov(w_i,w_j)为w_i与w_j的协方差；σ_wi与σ_wj为w_i与w_j的标准差；

建立满足给定相关性水平的输入变量概率模型的步骤为：

步骤二B：产生满足独立正态分布的样本矩阵S，其维数为m×N，N为样本规模；

步骤二C：运用三阶多项式正态变换(third-order polynomial normaltransformation，TPNT)以及求解代数方程的方法求取标准正态分布随机变量Y的等效相关系数矩阵R_Y；

步骤二D：对R_Y进行Cholesky分解获得其下三角矩阵L_Y，并由Y＝L_YS获得具有相关性的服从标准正态分布的样本矩阵Y；

步骤二E：通过等概率转换原则得到满足相关性水平的样本矩阵W，从而获得满足给定相关性水平的输入变量概率模型。

3.根据权利要求2所述的孤岛微电网概率潮流解析计算方法，其特征在于：步骤三中，所述基于孤岛微电网潮流模型进行基准运行点的稳态潮流计算包括以下步骤：

步骤三A：设定节点功率方程统一表达式：仅设定PQ节点与PV节点两种节点类型，其中PQ节点需列写有功与无功平衡方程，PV节点仅需列写有功平衡方程，方程形式为：

式中：P_Ci和Q_Ci分别为节点i恒功率电源注入的有功和无功；P_Di和Q_Di分别为节点i具有下垂控制特性的设备注入的有功和无功；P_Li和Q_Li分别为节点i的有功与无功负荷；P_i和Q_i分别为节点i注入的总有功功率和无功功率；

其中，所述具有下垂控制特性的设备注入的功率统一表示为：

(4)、(5)、(6)、(7)式中：P_Dimax、P_Dimin、Q_Dimax、Q_Dimin分别为节点i具有下垂控制特性的设备注入的有功的上限、节点i具有下垂控制特性的设备注入的有功的下限、节点i具有下垂控制特性的设备注入的无功的上限和节点i具有下垂控制特性的设备注入的无功的下限，若无具有下垂控制特性的设备则均取为0；f_max、f_min、U_max、U_min分别为系统频率允许的上限、系统频率允许的下限、电压允许的上限、电压允许的下限；K_Dfi、K_DUi分别为节点i对应的P-f、Q-U下垂系数；f₀、f、U₀、U_i分别为系统频率的空载值、系统频率、电压的空载值与电压在节点i的实际值；

其中，对负荷进行建模表示为：

(8)式中：P_LNi和Q_LNi分别为节点i在额定工况下的有功与无功负荷；U_Ni和f_N分别为节点i的额定电压与频率，取U_Ni＝1p.u.，f_N＝1p.u.；A_pi、B_pi、C_pi、A_qi、B_qi、C_qi分别为节点i的负荷有功功率中恒阻抗型的百分比系数、节点i的负荷有功功率中恒电流型的百分比系数、节点i的负荷有功功率中恒功率型的百分比系数、节点i的负荷无功功率中恒阻抗型的百分比系数、节点i的负荷无功功率中恒电流型的百分比系数、节点i的负荷无功功率中恒功率型的百分比系数，各满足总和为1；k_Lpi、k_Lqi分别为节点i负荷的有功和无功功率的静态频率调节系数；U_i为节点i系统电压的实际值；f为节点i系统频率的实际值；

其中，节点注入的有功与无功功率表示为：

(9)式中：n为节点数目；G_ij与B_ij分别为节点导纳矩阵的实部与虚部；δ_ij为节点i和节点j的相角差；U_i为节点i系统电压的实际值；U_j为节点j系统电压的实际值；

步骤三B：牛顿拉夫逊法求解:孤岛微电网中PQ节点与PV节点的数目分别为n_PQ与n_PV，则总计列写2n_PQ+n_PV个方程，采用牛顿拉夫逊法求解潮流方程组，其中，修正方程简写为：

(10)式中：ΔP、ΔQ为节点有功与无功不平衡量；Δf、Δδ、ΔU为分别为频率、相角、电压的修正量；J为雅克比矩阵，其分块矩阵分别为：

雅克比矩阵J中各分块矩阵的维数为：

各分块矩阵的表达式为：

(13)式中，各参数由(3)-(9)式中各参数代入计算。

4.根据权利要求3所述的孤岛微电网概率潮流解析计算方法，其特征在于：步骤四和步骤五的计算步骤具体为：

步骤四五A：在基准运行点处将潮流方程线性化并忽略2次及以上高次项，可得随机分量的关系为：

(14)式中：ΔW、ΔX、ΔZ分别为节点注入功率、节点状态变量以及支路潮流变量的随机分量；J₀为基准点处的雅克比矩阵；X₀为基准点处的节点状态量；

当已知W的各阶半不变量，则根据半不变量的线性性质获得X、Z的各阶半不变量；

步骤四五B：通过对R_W进行Cholesky分解获得其下三角矩阵L_W，将m个相关输入随机变量W转换为m个相互独立的变量W_ind的组合，转换关系如下：

W_ind＝L_W ^-1W (15)

相应的J₀中与原始m个相关输入随机变量W对应列的元素需要通过L_W修正；

步骤四五C：具有下垂控制特性的电源输出功率与系统频率及节点电压有关，为获取其输出的概率分布，需要计算具有下垂控制特性的设备输出功率的半不变量，假定在系统运行范围内，对式(4)(6)在基准点线性化，推得其关系为：

结合式(14)进一步转换为与节点注入功率ΔW的关系后同样基于半不变量的线性性质由W的各阶半不变量求得具有下垂控制特性的电源输出功率的各阶半不变量。

5.根据权利要求4所述的孤岛微电网概率潮流解析计算方法，其特征在于：步骤六中获取输出变量的概率密度函数的步骤为：

步骤六A：借助标准正态分布概率密度函数φ以及它的各阶导数的线性组合来逼近最终待求的分布函数，根据已求得状态变量的各阶半不变量为γ_k，列出由各阶半不变量推导各阶中心矩的公式：

步骤六B：由半不变量γ_k求得中心矩后，进一步通过中心矩表示Gram-Charlier级数系数，公式为：

步骤六C：将随机变量标准化处理后借助Gram-Charlier级数展开将随机变量的分布函数表达为由正态随机变量各阶导数组成的级数；

则由Gram-Charlier级数展开获取输出变量的概率密度函数及累积分布函数，输出概率潮流计算结果。

Images