一、代码整体概述
本套代码针对含风光发电的电力系统，考虑负荷波动及风、光出力不确定性，采用蒙特卡洛法与半不变量法（结合 Gram-Charlier 和 Cornish-Fisher 级数）实现概率潮流计算。代码基于 IEEE 34 节点和 IEEE 30 节点标准测试系统搭建，注释详尽，通过牛顿 - 拉夫逊法保证潮流计算收敛性，可输出节点电压、支路功率的概率分布及电压越限概率，为电力系统稳定性分析提供数据支撑。
二、核心功能模块
（一）基础数据与潮流计算模块
- 数据读取与初始化（dataIn.m、readdata.m）
- 读取 IEEE 34 节点（如 IEEE34_solar.txt）和 IEEE 30 节点（如 IEEE30.txt）系统参数，包括节点数、线路参数（电阻、电抗、电纳）、变压器参数（变比、阻抗）、电源与负荷参数（发电机出力、负荷有功/无功均值及标准差）。
- 初始化基准容量、最大迭代次数、计算精度等基础参数，区分 PQ 节点、PV 节点和平衡节点。
- 节点导纳矩阵构建（formACY.m）
- 基于线路和变压器参数，构建交流系统节点导纳矩阵 Y 及电纳矩阵 Y0，考虑线路充电电纳和变压器非标准变比的影响，采用稀疏矩阵存储以提升计算效率。
- 确定性潮流计算（NR_main.m、pf_line.m）
- 采用牛顿 - 拉夫逊法（NR）求解潮流，输出各节点电压幅值与相角、支路功率（有功/无功）及系统网损。
- 通过 formDetaPQ.m 计算有功、无功不平衡量，formJacco.m 生成雅可比矩阵，calcuVD.m 求解电压和相角修正量，确保迭代收敛。
（二）随机源建模模块
- 负荷波动建模
- 负荷服从正态分布，通过 textread 读取负荷数据（如 IEEE34load_30%.txt），获取各节点负荷有功/无功均值及标准差，在蒙特卡洛模拟中用 normrnd 函数生成随机负荷值。
- 光伏出力建模
- 基于广州或上海的光照强度数据（如广州光照强度数据样本.txt），采用 Beta 分布拟合光照强度概率特性，通过 HOMER 软件获取光强期望与方差，计算 Beta 分布参数 a、b。
- 结合光伏组件面积、光电转换效率，计算光伏随机出力的八阶原点矩与半不变量（NcalPLCum.m），并映射到对应接入节点（如 IEEE 34 节点的 34 号节点）。
- 风力发电建模
- 风速服从 Weibull 分布，通过 wblrnd 函数生成随机风速，结合风机特性（切入风速 vci、额定风速 vr、切出风速 vco），计算风力发电出力，接入指定节点（如 IEEE 30 节点的 15 号节点）。
（三）概率潮流计算模块
- 蒙特卡洛法（MC）实现（main_MonteCarlo.m、test.m 等）
- 设置模拟次数（默认 6000 次），每次迭代中随机生成负荷、风、光出力值，调用 NR_main.m 进行潮流计算，存储节点电压和支路功率数据。
- 统计节点电压、支路有功/无功的概率分布，通过 tabulate 函数生成频数表，计算累积概率，并采用 ksdensity 函数进行概率密度估计。
- 基于欧洲标准 EN 50160，通过 ProbMC.m 计算节点电压越限概率（电压＜0.95p.u.或＞1.05p.u.的概率）。
- 半不变量法实现（CM_MC_PLF.m、main_Gram_Charlier.m 等）
- 半不变量计算：通过 NcalGCum.m 计算发电机八阶半不变量，NcalPLCum.m 计算负荷和光伏的八阶半不变量，整合得到系统注入功率的半不变量矩阵 B_s。
- 灵敏度矩阵求解：基于雅可比矩阵的逆矩阵，计算电压和支路功率对注入功率的灵敏度矩阵（S1-S8、T1-T8）。
- 级数展开：
- Gram-Charlier 级数：通过 cv_s、cp_s、cq_s 计算电压、支路有功/无功的展开系数，结合正态分布概率密度函数，修正得到实际概率分布（f1_s、fp1_s、fq1_s）。
- Cornish-Fisher 级数：基于半不变量修正标准正态分布分位数，得到节点电压和支路功率的分位数分布（yyV、yyP、yyQ），计算电压越限概率（ProbCMCF.m）。
