基于蒙特卡洛与半不变量法的风光电力系统概率潮流计算算法研究——负荷波动与风光出力不确定性考量

含风光发电的电力系统概率潮流计算，考虑负荷波动，风力和光伏出力不确定性，算法方面：基于蒙特卡洛法和半不变量法（gram-charlier和corn-fisher级数）。 代码注释详尽，收敛性好。

一、代码整体概述

本套代码针对含风光发电的电力系统，考虑负荷波动及风、光出力不确定性，采用蒙特卡洛法与半不变量法（结合Gram-Charlier和Cornish-Fisher级数）实现概率潮流计算。代码基于IEEE 34节点和IEEE 30节点标准测试系统搭建，注释详尽，通过牛顿-拉夫逊法保证潮流计算收敛性，可输出节点电压、支路功率的概率分布及电压越限概率，为电力系统稳定性分析提供数据支撑。

二、核心功能模块

（一）基础数据与潮流计算模块

1. 数据读取与初始化（dataIn.m、readdata.m）
   - 读取IEEE 34节点（如IEEE34_solar.txt）和IEEE 30节点（如IEEE30.txt）系统参数，包括节点数、线路参数（电阻、电抗、电纳）、变压器参数（变比、阻抗）、电源与负荷参数（发电机出力、负荷有功/无功均值及标准差）。
   - 初始化基准容量、最大迭代次数、计算精度等基础参数，区分PQ节点、PV节点和平衡节点。
2. 节点导纳矩阵构建（formACY.m）
   - 基于线路和变压器参数，构建交流系统节点导纳矩阵Y及电纳矩阵Y0，考虑线路充电电纳和变压器非标准变比的影响，采用稀疏矩阵存储以提升计算效率。
3. 确定性潮流计算（NRmain.m、pfline.m）
   - 采用牛顿-拉夫逊法（NR）求解潮流，输出各节点电压幅值与相角、支路功率（有功/无功）及系统网损。
   - 通过formDetaPQ.m计算有功、无功不平衡量，formJacco.m生成雅可比矩阵，calcuVD.m求解电压和相角修正量，确保迭代收敛。

（二）随机源建模模块

1. 负荷波动建模
   - 负荷服从正态分布，通过textread读取负荷数据（如IEEE34load_30%.txt），获取各节点负荷有功/无功均值及标准差，在蒙特卡洛模拟中用normrnd函数生成随机负荷值。
2. 光伏出力建模
   - 基于广州或上海的光照强度数据（如广州光照强度数据样本.txt），采用Beta分布拟合光照强度概率特性，通过HOMER软件获取光强期望与方差，计算Beta分布参数a、b。
   - 结合光伏组件面积、光电转换效率，计算光伏随机出力的八阶原点矩与半不变量（NcalPLCum.m），并映射到对应接入节点（如IEEE 34节点的34号节点）。
3. 风力发电建模
   - 风速服从Weibull分布，通过wblrnd函数生成随机风速，结合风机特性（切入风速vci、额定风速vr、切出风速vco），计算风力发电出力，接入指定节点（如IEEE 30节点的15号节点）。

（三）概率潮流计算模块

1. 蒙特卡洛法（MC）实现（mainMonteCarlo.m、test.m等）
   - 设置模拟次数（默认6000次），每次迭代中随机生成负荷、风、光出力值，调用NRmain.m进行潮流计算，存储节点电压和支路功率数据。
   - 统计节点电压、支路有功/无功的概率分布，通过tabulate函数生成频数表，计算累积概率，并采用ksdensity函数进行概率密度估计。
   - 基于欧洲标准EN 50160，通过ProbMC.m计算节点电压越限概率（电压＜0.95p.u.或＞1.05p.u.的概率）。
2. 半不变量法实现（CMMCPLF.m、mainGramCharlier.m等）
   -半不变量计算：通过NcalGCum.m计算发电机八阶半不变量，NcalPLCum.m计算负荷和光伏的八阶半不变量，整合得到系统注入功率的半不变量矩阵Bs。
   -灵敏度矩阵求解：基于雅可比矩阵的逆矩阵，计算电压和支路功率对注入功率的灵敏度矩阵（S1-S8、T1-T8）。
   -级数展开：
   - Gram-Charlier级数：通过cvs、cps、cqs计算电压、支路有功/无功的展开系数，结合正态分布概率密度函数，修正得到实际概率分布（f1s、fp1s、fq1_s）。
   - Cornish-Fisher级数：基于半不变量修正标准正态分布分位数，得到节点电压和支路功率的分位数分布（yyV、yyP、yyQ），计算电压越限概率（ProbCMCF.m）。

（四）结果分析与可视化模块

1. 结果对比：通过difference1、difference2计算半不变量法与蒙特卡洛法的电压越限概率误差，验证半不变量法的准确性。
2. 图形绘制：
   - 概率密度曲线：绘制节点电压（如节点33、34）、支路功率（如线路31-33）的概率密度曲线，对比蒙特卡洛法直方图、半不变量法解析结果。
   - 累积分布曲线：展示不同方法下节点电压和支路功率的累积概率，直观呈现随机变量分布特性。
   - 越限概率统计：输出各PQ节点的电压越限概率，判断系统运行风险。

三、关键技术特点

1. 不确定性建模全面：涵盖负荷波动、光伏（Beta分布）和风电（Weibull分布）出力不确定性，贴合实际电力系统特性。
2. 算法高效可靠：蒙特卡洛法保证计算精度，半不变量法结合级数展开大幅提升计算速度，适用于大规模系统。
3. 收敛性保障：牛顿-拉夫逊法迭代求解潮流，通过设置最大迭代次数和精度阈值，确保潮流计算收敛；半不变量法通过八阶展开减少截断误差。
4. 通用性强：支持IEEE 34/30节点系统，可通过修改数据文件适配不同电网结构，灵活调整风、光接入节点及容量。

四、使用说明

1. 数据准备：确保IEEE34solar.txt、IEEE30.txt等系统参数文件，以及负荷、光照强度数据文件（如IEEE34load30%.txt、广州光照强度数据样本.txt）路径正确。
2. 参数调整：根据需求修改模拟次数（daishu）、光伏组件面积（A）、光电转换效率（prey）、风机特性参数（vci、vr、vco）等。
3. 运行流程：
   - 先运行dataIn.m读取系统数据，再通过NRmain.m验证确定性潮流收敛性。
   - 选择算法：运行mainMonteCarlo.m执行蒙特卡洛法，运行CMMCPLF.m执行半不变量法（含两种级数展开）。
   - 查看结果：通过图形窗口查看概率密度和累积分布曲线，查看命令行输出的电压越限概率和计算时间。

五、典型应用场景

1. 含高比例新能源的配电网规划：分析风光接入对节点电压稳定性的影响，优化风光装机容量和接入位置。
2. 电力系统风险评估：计算不同运行场景下的电压越限概率，为系统调度提供决策依据。
3. 算法对比研究：作为概率潮流算法验证的基准代码，对比蒙特卡洛法与半不变量法的精度和效率。

要不要我帮你整理一份代码核心函数调用关系图？可以直观展示各模块间的数据流和函数依赖，方便快速理解代码逻辑。