【计及DG和相关性的纯交流三点估计随机潮流计算方法】

目标需求##

计及DG和相关性的纯交流三点估计随机潮流(matlab版)

目标问题

1. 要适用于任意大小的纯交流电网,支持节点和支路的增删；
2. 要考虑了DG,负荷波动情况，并能够处理相关性计算
3. 考虑采用Gram_Charlier级数或CornishFish级数拟合概率分布函数；
4. 如何构建原始状态变量X与标准独立状态空间的变量Z转换；
5. 如何验证方法误差在有效的范围内，考虑采用蒙特卡洛作对比；

初定子函数列写清单

1. 主程序命名：main2_2
2. 输入算例命名：data_ieee33
3. 牛拉法基波子命名：runpf()(提前安装matpower)
4. 节点导纳子命名：createYbus
5. 雅克比矩阵子命名：Jacobi(n-1+m*n-1+m)
6. 相关性Nataf逆变换命名：NatafMethod
7. 变量X与Y之间的转换：EqualYX文件(查阅参考文献)
   
   建立类似的工作文档，一目了然
   

程序计算步骤及流程图

以matlab版程序作为讲解。

1 参数初始化

参数初始化，读取网络参数，在此标幺化，随机潮流是在确定性基波潮流计算基础之上进行的，参数初始化和基波潮流学习： https://blog.csdn.net/WConstelltion/article/details/123751611.

2 蒙特卡洛计算

2.1负荷抽样：(蒙特卡洛相关性抽样)

	sampleZ=mvnrnd(muz,sigmaz,m)'; 
	sampleY = B*sampleZ; %构建相关正态随机变量
	[sampleX,muX,sigmaX] = EqualYX(Prob,sampleY);

2.2DG抽样：(风速、光强、储能功率)

     v_FD = wblrnd(FD_c,FD_k,m,1); %风速抽样
    cd_GF= betarnd(GF_a,GF_b,1,m);   %光强抽样

2.3相关随机变量样本矩阵转为节点功率矩阵：

mc_point_vbQ(i, :)  = load_tan(i).*mc_point_vbP(i, :);%负荷点
mc_point_vbQ(n_load+i,:) = tan(acos(lanta))*mc_point_vbP(n_load+i,:);%风电
mc_point_vbP(n_load+i,:) = mc_point_vbP(n_load+i,:)*yita*A/1000;%光伏

2.4节点等效负荷矩阵：

    vb_dgi = sum(index_load(1:dg(i,2)));%第i个DG接入节点所对应变量矩阵的行号
    mc_point_vbP(vb_dgi, :) = mc_point_vbP(vb_dgi, :)-mc_point_vbP(n_load+i, :);
    mc_point_vbQ(vb_dgi, :) = mc_point_vbQ(vb_dgi, :)-mc_point_vbQ(n_load+i, :);%节点综合注入功率

2.5修改节点负荷参数：

        shuju.bus(index_load,3)=load_p(:,i);
        shuju.bus(index_load,4)=load_q(:,i);

2.6确定性潮流计算：

       [basemva,bus,gen,branch]=runpf(shuju, mpopt);

3 基于Nataf逆变换的3PEM计算流程

本文算法流程引用前辈的计算思路并基于改进：https://download.csdn.net/download/destiny9613/10004292?spm=1001.2101.3001.6650.17&utm_medium=distribute.pc_relevant.none-task-download-2%7Edefault%7EBlogCommendFromBaidu%7ERate-17.pc_relevant_paycolumn_v3&depth_1-utm_source=distribute.pc_relevant.none-task-download-2%7Edefault%7EBlogCommendFromBaidu%7ERate-17.pc_relevant_paycolumn_v3&utm_relevant_index=24.
3.1求解随机变量相关系数矩阵的下三角矩阵B：

[rhoY,B] = NatafMethod(Prob);%相关系数

3.2随机变量独立空间样本点矩阵Z：

location_zp=[muZ muZ muZ]+[zeta1_p zeta2_p zeta3_p].*[sigmaZ sigmaZ sigmaZ];

3.3考虑相关性以及等概率变换过程

location_yp=B*location_zp; %去XY相关性
[location_xp,muX,sigmaX] = EqualYX(Prob,location_yp);

3.4样本点的权重计算

weight1=1./(zeta1_p.*(zeta1_p-zeta2_p));
weight2=-1./(zeta2_p.*(zeta1_p-zeta2_p));
weight3=1/Prob.Nx-1./(lambda4_p-lambda3_p.^2);

3.5变量估计点转化为对应的节点功率点矩阵(包含了DG功率)

point_vbP = repmat(location_p(:, 3),1,2*Prob.Nx+1);%形成nx*2nx+1维随机变量期望矩阵
    point_vbP(i,2*i-1)=location_p(i,1);%nx*3为矩阵
    point_vbP(i,2*i)=location_p(i,2);

3.6形成节点功率估计点矩阵

point_p=point_vbP(1:n_load, :); %节点负荷赋值
point_q = point_vbQ(1:n_load, :); 
point=[point_p;point_q]; 

3.7蒙特卡洛节点电压概率分布

mc_vm_mu=mean(mc_vm,2);%返回矩阵各行均值
mc_vm_pdf(i, :) = ksdensity(mc_vm(i, :),vm_xi,'function','pdf'); %离散点求概率密度函数

3.8CornishFish展开系数
以节点电压为例分析

pem_vm_mu= aerfa_vm(:, 1);    %一阶原点矩
pem_vm_beta= aerfa_vm(:, 2);  %二阶原点矩
pem_vm_std=sqrt(pem_vm_beta-pem_vm_mu.^2);  
[pem_vm_pdf,pem_vm_cdf]= CornishFish(aerfa_vm,range_vm,pem_vm_mu,pem_vm_std);  %级数展开求解概率分布

节点电压计算结果误差相对较小

程序连接

铁子们点它: https://download.csdn.net/download/lisha1002/12148929?ops_request_misc=%257B%2522request%255Fid%2522%253A%2522164835259216780366543218%2522%252C%2522scm%2522%253A%252220140713.130102334.pc%255Fdownload.%2522%257D&request_id=164835259216780366543218&biz_id=1&utm_medium=distribute.pc_search_result.none-task-download-2_downloadfirst_rank_ecpm_v1~rank_v31_ecpm-5-12148929.pc_v2_rank_dl_default&utm_term=%E5%8D%8A%E4%B8%8D%E5%8F%98%E9%87%8F&spm=1018.2226.3001.4451.5.
编写不易，需要考虑多种因素的程序私信我铁子们(网站拿来主义多)……