目录
主要内容
该程序主要内容是基于半不变量法的概率潮流，包含蒙特卡洛模拟法、半不变量法＋Gram-Charlier级数展开以及半不变量法+Cornish-Fisher级数展开三种方法以及效果对比，模型考虑了负荷不确定性（正态分布）和光伏随机特性，以IEEE34节点为例，计算节点电压、支路潮流概率密度、累计概率并绘制曲线。程序采用matlab编写，注释清晰，主函数为main_CM_MC，子程序较多，适合有经验的人入手学习！
部分代码
%===================考虑分布式电源、发电机和负荷随机波动的概率潮流计算================================
tic %半不变量法计算计时开始
%% 基础参数------------------------------------------------------------------
[Nodes,linenum,SB,maxIters,OPdata1,precision,OPdata2,balanceID,balancenotes,...
lineID,linei,linej,liner,linex,lineb,...
branchi,branchb,...
transID,transi,transj,transr,transx,transk,transkMin,transkMax,...
PQi,PG,QG,PD,QD,...
PVi,PVV,PVQmin,PVQmax...
NGi,OP_0,OP_1,OP_2,NGmin,NGmax]=dataIn('IEEE34.txt'); %% 将数据放入各变量后以列向量的格式输出
%% 首先进行基础潮流计算，形成雅克比矩阵
%形成交流系统节点导纳矩阵----------------------------------------------------
[Y,Y0] = formACY(Nodes,branchi,branchb,linei,linej,liner,...
linex,lineb,transi,transj,transr,transx,transk);
%潮流计算-------------------------------------------------------------------
[V,deta,PQ_loss,S,detaS,Colab,Jacco,Jacco2 ]...
= NR_main(PVi,PVV,balancenotes,Y,Y0,linei,linej,transi,transj,...
PG,PD,QG,QD,maxIters,precision,Nodes);
%% 计算输入的半不变量
%发电机的随机参数输入--------------------------------------------------------
%%pdfgen(i,1)为发电机序号
%%pdfgen(i,2)为发电机的节点号
%%pdfgen(i,3)为发电机的有功出力
%%pdfgen(i,4)为发电机的无功出力
%%pdfgen(i,5)为发电机的出力的概率
%发电机的八阶半不变量形成-------------------------------------------------------
pdfgen=textread('IEEE34gen.txt');%%普通发电机出力服从二项分布;
ngen=length(pdfgen(:,1));
PgPx=zeros(Nodes,8);
PgQx=zeros(Nodes,8);
PgPx(pdfgen(:,2),:)=NcalGCum(pdfgen(:,3),pdfgen(:,5));
PgQx(pdfgen(:,2),:)=NcalGCum(pdfgen(:,4),pdfgen(:,5));
%负荷的八阶半不变量形成-------------------------------------------------------
%%pdfload(i,1)为负荷序号
%%pdfload(i,2)为负荷的节点号
%%pdfload(i,3)为负荷有功均值
%%pdfload(i,4)为负荷无功均值
%%pdfload(i,5)为负荷有功标准差
%%pdfload(i,6)为负荷无功标准差 %%标准差给定可以参照“3Sita原则”
%负荷的八阶半不变量---------------------------------------------------------
pdfload=textread('IEEE34load_30%.txt');%%负荷负荷正态分布
nload=length(pdfload(:,1));
PlPx=zeros(Nodes,8);
PlQx=zeros(Nodes,8);
PlPx(pdfload(:,2),:)=NcalPLCum(-pdfload(:,3),-pdfload(:,5));
PlQx(pdfload(:,2),:)=NcalPLCum(-pdfload(:,4),-pdfload(:,6));
%-------光伏随机特性建模-----------------------------
%选择上海31°8’N、121°35’E作为光照强度分布的考量位置，在HOMERE软件上获取光强分布的期望值和方差。
%miu=0.150314263;
%sita=0.049758487;
%利用HOMER软件获取广州（113°15′E，23°7′N）的光照强度数据样本作为后续应用的模型
M=textread('Guangzhao.txt');