MATLAB配电网分布式电源承载力评估框架

简介本资源是《新型电力系统下多分布式电源接入配电网承载力评估方法研究》论文的完整MATLAB复现代码包面向电力系统方向的研究生、科研人员及工程技术人员解决分布式电源如光伏、风电大规模接入后配电网承载能力量化评估难题。压缩包共11个文件含9个核心MATLAB脚本如main.m主程序、AHP.m与entropy_weight.m分别实现层次分析法与熵权法组合赋权、cal_index.m计算评价指标等、1份PDF代码说明文档和1个Excel风光负荷曲线数据表总大小仅1.31MB轻量易部署。已有1291人学习下载代码结构清晰、模块功能明确覆盖二阶锥规划建模、IEEE33节点系统仿真、多维度承载力指标计算及综合评分生成全流程配套说明文档详述各函数作用与调用逻辑便于快速理解算法原理、验证评估结果并开展二次开发。1. 分布式电源接入配电网承载力评估不是算个潮流就完事而是要回答“还能接多少、在哪接、怎么接才不越限”你手头有一片10kV辐射状配电网已有3台光伏逆变器在馈线中段运行最近又有两个屋顶光伏项目申请并网——调度中心甩来一句“评估下承载力一周内给结论。”你打开MATLAB翻出去年写的潮流计算脚本跑完一看电压偏差0.8%支路负载率72%一切正常。于是你信心满满回复“可接入”。结果现场调试时某节点夜间无功倒送导致电容器频繁投切白天光伏大发又引发末端电压越上限保护装置误动作两次。这不是算得不准是承载力评估被窄化成了单点稳态校核。真正的分布式电源承载力评估是面向规划与运行双场景的多维度约束综合判断既要守住电压偏差、短路电流、谐波畸变、反向功率、保护配合这五条硬边界又要量化不同接入位置、容量组合、控制策略下的裕度衰减曲线。它不依赖某套“万能公式”而是一套可配置、可追溯、可复用的MATLAB工程化分析框架——本文复现的正是这样一套落地路径从原始网架建模、多工况潮流扫描、约束敏感性排序到承载力热力图生成全部基于原生MATLAB无需Simulink或第三方工具箱适配2021b及以上版本代码结构清晰、注释完整、参数接口开放。如果你正被“评估报告被退回三次”“模型改一次就要重跑半天”“领导问‘为什么这个位置比那个位置多接50kW’答不上来”困扰这篇就是为你写的。2. 用MATLAB构建可复用的配电网承载力评估框架从网架建模到多工况扫描承载力评估不是一次性快照而是一场有组织的“压力测试”。核心在于建立一个可参数化、可迭代、可验证的MATLAB分析框架。我们不堆砌复杂算法而是用最扎实的原生功能搭起四层骨架网架数据层 → 潮流引擎层 → 约束校核层 → 承载力映射层。每层都设计为独立函数输入输出明确便于替换算法或扩展指标。2.1 网架建模用结构体替代Excel表格让拓扑关系可编程很多工程师习惯把线路参数、节点负荷、DG位置全塞进Excel每次换网架就手动改表——这是复现失败的第一大坑。正确做法是用MATLAB结构体struct定义网架将物理关系转化为代码逻辑% 定义配电网拓扑结构体示例IEEE 33节点简化版 net struct(); net.nbus 33; % 节点总数 net.nline 32; % 支路总数 net.baseMVA 10; % 基准容量MVA net.baseKV 12.66; % 基准电压kV % 节点数据[节点ID, 类型(1PV,2PQ,3Slack), P0(MW), Q0(MVar), Vmin(pu), Vmax(pu)] net.bus [ 1, 3, 0, 0, 0.95, 1.05; % 平衡节点 2, 2, 0.1, 0.06, 0.95, 1.05; % PQ负荷节点 3, 2, 0.09, 0.04, 0.95, 1.05; % ... 其余30行省略实际需填满33行 ]; % 支路数据[支路ID, 送端, 受端, R(pu), X(pu), B/2(pu), Imax(pu)] net.line [ 1, 1, 2, 0.0922, 0.047, 0, 1.2; 2, 2, 3, 0.493, 0.2511, 0, 1.2; % ... 其余31行 ]; % 分布式电源接入点预设支持多点、多类型 net.dg_sites [2, 5, 10, 15, 20]; % 可接入节点ID列表 net.dg_types {PV, WT, BESS}; % 支持的DG类型逻辑说明net结构体封装了全部拓扑与电气参数net.bus中第5、6列直接定义电压约束上下限net.line中最后一列Imax即支路热稳定极限标幺值。这种写法让网架成为“活数据”——后续修改只需调整结构体字段无需碰任何计算逻辑。参数说明baseMVA和baseKV决定标幺化基准必须与实际系统一致Vmin/Vmax取值参考DL/T 1648-2016《配电网接纳分布式电源能力评估导则》通常设为0.95–1.05 puImax按电缆长期允许载流量折算例如YJV22-3×240mm²在空气中敷设约400A对应标幺值需按baseI baseMVA / (sqrt(3)*baseKV)换算。2.2 潮流引擎手写前推回代法比调用loadflow更可控Simulink或Power System Toolbox的潮流函数虽快但内部黑匣子无法干预收敛过程、无法注入DG控制逻辑、无法获取每步迭代的中间变量。承载力评估恰恰需要这些——比如观察“当光伏出力从0升至额定值时某节点电压如何连续变化”。因此我们采用手写前推回代Backward/Forward Sweep潮流算法代码仅120行却完全透明function [V, Sloss] backward_forward_sweep(net, Pdg, Qdg) % 输入net-网架结构体Pdg/Qdg-各节点DG有功/无功出力向量长度net.nbus % 输出V-节点电压复数向量Sloss-总损耗复数 V ones(net.nbus, 1) * (1 1j*0); % 初始电压全设为1∠0° max_iter 20; tol 1e-5; for iter 1:max_iter % --- 回代从末端向首端计算支路电流 --- Iline zeros(net.nline, 1); % 支路电流复数向量 for i net.nline:-1:1 from net.line(i,2); to net.line(i,3); % 节点注入功率 负荷 - DG出力 下游支路流入功率 S_node (net.bus(to,3) - Pdg(to)) 1j*(net.bus(to,4) - Qdg(to)); if to net.nbus % 非末端节点加上下游支路功率 idx_down find(net.line(:,2) to); if ~isempty(idx_down), S_node S_node sum(Iline(idx_down).*conj(V(to))); end end Iline(i) conj(S_node / V(to)); % 该支路电流 注入功率 / 电压共轭 end % --- 前推从首端向末端更新节点电压 --- V_new V; for i 1:net.nline from net.line(i,2); to net.line(i,3); Z_line net.line(i,4) 1j*net.line(i,5); % 支路阻抗 V_new(to) V_new(from) - Iline(i) * Z_line; end % --- 收敛判断 --- if norm(abs(V_new) - abs(V), inf) tol V V_new; break; end V V_new; end % 计算总损耗 Sloss 0; for i 1:net.nline from net.line(i,2); to net.line(i,3); I Iline(i); V_from V(from); Sloss Sloss I * conj(I) * (net.line(i,4) 1j*net.line(i,5)); end end逻辑说明该函数严格遵循配电网辐射状特性回代阶段从叶节点向上累加注入功率前推阶段从根节点向下逐级压降。关键创新点在于Iline计算中显式分离了“负荷消耗”与“DG注入”为后续承载力扫描预留了接口——Pdg和Qdg向量可动态赋值实现任意DG组合的快速遍历。参数说明Pdg/Qdg必须是长度为net.nbus的向量未接入DG的节点对应位置填0tol1e-5确保电压幅值误差小于0.001pu对承载力评估足够若遇不收敛优先检查net.line中支路方向是否全部为“送端→受端”即from to这是前推回代法的前提。2.3 多工况扫描用嵌套for循环生成承载力搜索空间拒绝暴力穷举承载力不是“最大能接多少kW”而是“在满足所有约束前提下各接入点的最大可行容量组合”。暴力穷举如对5个点各试100个容量值需10^10次计算不可行。我们采用分层扫描策略先固定DG类型与位置对单点做容量扫描得到“单点承载力曲线”再基于曲线斜率筛选高敏感节点对其做联合扫描。核心是scan_capacity函数function [cap_curve, constraint_violation] scan_capacity(net, dg_site, dg_type, cap_step, cap_max) % 输入dg_site-接入节点IDdg_type-PV/WT/BESScap_step-步长(kW)cap_max-上限(kW) % 输出cap_curve-容量-约束裕度矩阵constraint_violation-各约束越限标志矩阵 n_steps floor(cap_max / cap_step) 1; cap_curve zeros(n_steps, 6); % 列[容量, Vmin, Vmax, Imax, Ishort, THD] constraint_violation false(n_steps, 5); % 行各约束是否越限 for k 1:n_steps cap_kW (k-1) * cap_step; Pdg zeros(net.nbus, 1); Qdg zeros(net.nbus, 1); % 根据DG类型设置出力模型简化版实际需接气象/负荷数据 switch dg_type case PV Pdg(dg_site) cap_kW / 1000; % 转为MW Qdg(dg_site) 0; % 恒功率因数1 case WT Pdg(dg_site) cap_kW / 1000 * 0.3; % 风电容量系数取0.3 Qdg(dg_site) 0; case BESS Pdg(dg_site) cap_kW / 1000; % 储能充放电功率 Qdg(dg_site) 0; end % 调用潮流计算 [V, ~] backward_forward_sweep(net, Pdg, Qdg); % 校核五大约束简化版实际需调用详细函数 cap_curve(k,1) cap_kW; cap_curve(k,2) min(abs(V)); % 最小电压幅值(pu) cap_curve(k,3) max(abs(V)); % 最大电压幅值(pu) cap_curve(k,4) max_line_loading(net, V, Pdg, Qdg); % 最大支路负载率 cap_curve(k,5) short_circuit_current(net, dg_site, cap_kW); % 短路电流(kA) cap_curve(k,6) harmonic_distortion(net, dg_site, cap_kW); % 谐波畸变率(%) % 标记越限 constraint_violation(k,1) cap_curve(k,2) net.bus(dg_site,5); % 电压下限 constraint_violation(k,2) cap_curve(k,3) net.bus(dg_site,6); % 电压上限 constraint_violation(k,3) cap_curve(k,4) 1.0; % 支路过载 constraint_violation(k,4) cap_curve(k,5) 20; % 短路电流限值 constraint_violation(k,5) cap_curve(k,6) 5; % 谐波限值(5%) end end % 辅助函数计算最大支路负载率 function Imax_pu max_line_loading(net, V, Pdg, Qdg) Imax_pu 0; for i 1:net.nline from net.line(i,2); to net.line(i,3); Z net.line(i,4) 1j*net.line(i,5); I (Pdg(to) - 1j*Qdg(to)) / conj(V(to)); % 简化忽略线路损耗 Imax_pu max(Imax_pu, abs(I) / net.line(i,6)); end end逻辑说明scan_capacity对单点生成一条“容量-约束响应曲线”每行记录该容量下5个关键约束的实时值。constraint_violation矩阵是布尔型直接标出哪一容量步触发了哪个约束——这才是承载力的真正答案承载力由最先越限的约束决定而非某个平均值。例如某曲线显示容量达120kW时电压越上限约束2为true而支路负载率仅0.85此时承载力就是120kW且瓶颈是电压调节能力。参数说明cap_step建议取5–10kW过大会漏掉拐点cap_max设为预估上限的1.5倍如规划容量200kW则设300kWdg_type影响出力模型PV按日照模型、WT按风速模型、BESS需考虑SOC约束此处为简化用恒定系数实际项目中应替换为pv_power_model.m等专用函数。3. 承载力瓶颈识别与可视化用热力图定位“脆弱节点”告别模糊结论评估报告里写“承载力不足”毫无价值决策者需要知道“哪个约束在哪个节点最先卡住脖子” 这要求我们超越单点扫描构建全网承载力空间映射。核心是grid_capacity_map函数它对每个可接入节点执行scan_capacity提取其“临界容量”首个越限点对应的容量再插值生成二维热力图。3.1 生成全网承载力热力图地理坐标驱动的插值映射配电网承载力强弱与地理布局强相关——主干线上游节点通常裕度大末端分支节点极易越限。因此热力图必须基于真实地理坐标而非节点编号顺序function capacity_map grid_capacity_map(net, dg_types, cap_step, cap_max) % 输入net-含地理坐标的网架结构体dg_types-DG类型列表 % 输出capacity_map-地理网格上的承载力矩阵单位kW % 步骤1提取所有可接入节点的地理坐标需提前录入net.bus_geo % 假设net.bus_geo为nbus×2矩阵[x,y]单位米 x_coords net.bus_geo(:,1); y_coords net.bus_geo(:,2); % 步骤2对每个dg_site执行scan_capacity获取临界容量 crit_cap zeros(net.nbus, 1); for i 1:net.nbus if ismember(i, net.dg_sites) % 仅对允许接入的节点计算 [~, viol] scan_capacity(net, i, dg_types{1}, cap_step, cap_max); % 找到第一个越限的容量步 first_viol find(any(viol,2), 1, first); if ~isempty(first_viol) crit_cap(i) (first_viol-1) * cap_step; else crit_cap(i) cap_max; % 全程未越限 end end end % 步骤3地理插值生成热力图使用MATLAB内置scatteredInterpolant F scatteredInterpolant(x_coords, y_coords, crit_cap, natural); % 定义绘图区域网格覆盖整个配电网地理范围 xq linspace(min(x_coords), max(x_coords), 200); yq linspace(min(y_coords), max(y_coords), 200); [Xq, Yq] meshgrid(xq, yq); capacity_map F(Xq, Yq); % 步骤4掩膜处理——只显示配电网覆盖区域需提供地理边界poly if isfield(net, boundary_poly) in_boundary inpolygon(Xq(:), Yq(:), net.boundary_poly(:,1), net.boundary_poly(:,2)); capacity_map(~reshape(in_boundary, size(Xq))) NaN; end end逻辑说明scatteredInterpolant采用自然邻域插值natural比线性插值更能反映配电网的拓扑突变特性——例如两条平行馈线间承载力可能骤降50%线性插值会平滑掉这种关键梯度。capacity_map矩阵中每个像素值代表该地理位置的“等效承载力”单位kW直接对应规划人员最关心的“此处最多能装多少”。参数说明net.bus_geo必须在建模阶段录入可通过GIS软件导出CSV后用readmatrix加载net.boundary_poly是配电网服务区域的多边形顶点坐标用于裁剪热力图避免空白区域干扰判断插值网格分辨率200×200已足够更高分辨率仅增加计算量不提升决策精度。3.2 瓶颈约束溯源用堆叠柱状图揭示“谁在拖后腿”热力图告诉我们在哪弱堆叠柱状图则解释为什么弱。对任一节点我们绘制其5个约束的“裕度衰减曲线”直观显示哪个约束最先触顶function plot_constraint_bottleneck(net, dg_site, dg_type, cap_step, cap_max) [cap_curve, viol] scan_capacity(net, dg_site, dg_type, cap_step, cap_max); % 计算各约束裕度正值为裕度负值为越限量 margin_Vmin cap_curve(:,2) - net.bus(dg_site,5); % 电压下限裕度 margin_Vmax net.bus(dg_site,6) - cap_curve(:,3); % 电压上限裕度 margin_Imax 1.0 - cap_curve(:,4); % 支路负载率裕度 margin_Ishort 20 - cap_curve(:,5); % 短路电流裕度(kA) margin_THD 5 - cap_curve(:,6); % 谐波裕度(%) % 绘制堆叠图横轴为容量纵轴为各约束裕度 figure(Position,[100,100,800,500]); stackplot(cap_curve(:,1), margin_Vmin, margin_Vmax, ... margin_Imax, margin_Ishort, margin_THD, ... DisplayName,{电压下限,电压上限,支路负载,短路电流,谐波畸变}); xlabel(分布式电源容量 (kW)); ylabel(约束裕度); title(sprintf(节点%d承载力瓶颈分析%s, dg_site, dg_type)); legend(Location,southoutside,Orientation,horizontal); grid on; % 标出首个越限点 first_viol_idx find(any(viol,2), 1, first); if ~isempty(first_viol_idx) hold on; plot(cap_curve(first_viol_idx,1), 0, ro, MarkerSize,10, LineWidth,2); text(cap_curve(first_viol_idx,1), 0.1, sprintf(临界点:%dkW, cap_curve(first_viol_idx,1)), ... HorizontalAlignment,center,FontSize,10,FontWeight,bold); end end逻辑说明stackplot将5个约束的裕度堆叠在同一坐标系曲线从上至下依次代表不同约束。当某条曲线首次跌穿横轴y0即表示该约束在此容量下越限。图中红色标记点直指承载力天花板且堆叠顺序天然暴露了“主次矛盾”——例如若电压上限曲线始终在最底层说明该节点电压调节能力是绝对瓶颈优化方向应聚焦于SVG或OLTC配置而非纠结于谐波滤波器。参数说明margin_*计算中net.bus(dg_site,5/6)直接读取节点定义的电压限值确保与建模一致stackplot是MATLAB R2022a新增函数若用旧版本可用barhold on手动堆叠图中southoutside图例位置避免遮挡关键曲线。3.3 承载力报告自动生成用MATLAB Report Generator输出PDF告别Word粘贴评估成果必须形成可交付、可追溯、可审计的正式报告。MATLAB Report Generator需安装能将分析结果自动编译为专业PDF包含图表、数据表、结论摘要function generate_capacity_report(net, capacity_map, dg_site_list, report_name) import mlreportgen.report.*; import mlreportgen.dom.*; rpt Report(report_name, pdf); add(rpt, TitlePage(Title,分布式电源接入配电网承载力评估报告,... Subtitle,[日期 datestr(now) , 网架 net.name])); % 执行摘要页 add(rpt, Chapter(执行摘要)); add(rpt, Paragraph(本报告基于MATLAB平台对 net.name 配电网开展分布式电源承载力评估。)); add(rpt, Paragraph([全网平均承载力 num2str(mean_capacity(capacity_map)) kW])); add(rpt, Paragraph([承载力最低节点 num2str(min_node_id) num2str(min_capacity) kW])); % 热力图页 add(rpt, Chapter(全网承载力空间分布)); img Image(capacity_heatmap.png); img.Width 6in; img.Height 4in; add(rpt, img); add(rpt, Paragraph(图承载力热力图单位kW。颜色越深表示承载力越高。)); % 关键节点详情页 for i 1:min(length(dg_site_list),3) % 仅展示前3个关键节点 add(rpt, Chapter([节点 num2str(dg_site_list(i)) 详细分析])); % 插入该节点的瓶颈堆叠图 fig_handle plot_constraint_bottleneck(net, dg_site_list(i), PV, 5, 300); saveas(fig_handle, [node_ num2str(dg_site_list(i)) _bottleneck.png]); close(fig_handle); img2 Image([node_ num2str(dg_site_list(i)) _bottleneck.png]); img2.Width 6in; img2.Height 4in; add(rpt, img2); end close(rpt); rptview(report_name .pdf); % 自动打开PDF end逻辑说明Report Generator将MATLAB计算结果与排版逻辑解耦——图表生成、数据提取、文字撰写全部由代码驱动杜绝人工粘贴导致的数值错位。add(rpt, TitlePage(...))自动生成带日期的封面Chapter自动编号Paragraph支持LaTeX数学公式如插入$P_{DG}^{max}120\text{kW}$。最终PDF可直接提交给调度中心符合DL/T 1648-2016对评估报告的格式要求。参数说明report_name建议包含网架名与日期如XX变电站10kV馈线_20240520.pdfmin_node_id/min_capacity需在调用前计算体现报告严谨性saveas保存图表时指定.png格式确保PDF嵌入清晰。4. 避坑指南承载力评估中5个血泪经验总结少踩一个省三天承载力评估看似流程化实则处处是坑。以下5条均来自真实项目翻车现场每一条都附带现象、原因与可立即执行的解决动作。请务必在启动新项目前逐条核对。4.1 现象潮流计算收敛但承载力曲线在低容量区就出现电压越限原因网架结构体net.bus中节点类型第2列设置错误。将本应为PQ负荷节点误设为PV节点导致潮流程序强制维持该节点电压为1.0pu掩盖了真实电压分布。解决在net.bus初始化后添加校验代码% 强制校验除平衡节点外所有DG接入点及负荷节点必须为PQ型类型2 pq_nodes find(net.bus(:,2) 2); slack_nodes find(net.bus(:,2) 3); if length(slack_nodes) ~ 1 || any(pq_nodes 1) % 节点1必须是平衡节点 error(网架节点类型定义错误仅节点1可为Slack其余必须为PQ); end4.2 现象同一节点PV与WT承载力评估结果差异巨大但现场实测两者表现接近原因scan_capacity中DG出力模型过于简化。PV按100%出力、WT按30%容量系数未考虑二者出力时序互补性——实际中光伏午间大发、风电夜间出力叠加后系统压力反而小于单一类型。解决改用时序联合扫描。准备典型日出力曲线如pv_profile_24h.csv、wt_profile_24h.csv在scan_capacity中改为% 加载24小时出力曲线每行1小时列数节点数 pv_prof readmatrix(pv_profile_24h.csv); % 33×24矩阵 wt_prof readmatrix(wt_profile_24h.csv); % 计算24小时中最大电压偏差、最大支路负载率 max_Vdev max(abs(V) - 1); % 取24小时最严重时刻4.3 现象热力图显示某区域承载力极高但现场该处变压器已重载原因网架建模遗漏了配电变压器。net.line只包含10kV线路未包含配变如S11-500kVA及其低压侧阻抗导致潮流计算未计入配变损耗与压降高估了末端节点电压裕度。解决将配变建模为“虚拟支路”。在net.line末尾追加% 在节点i后接入配变虚拟支路ID33 net.line(end1,:) [33, i, i1, 0.002, 0.015, 0, 1.2]; % R,X按500kVA折算 % 新增低压节点i1类型为PQ负荷接在其上 net.bus(end1,:) [i1, 2, 0.3, 0.15, 0.95, 1.05];4.4 现象短路电流计算值远低于实测值承载力评估过于乐观原因短路电流计算未计入DG提供的短路电流。传统方法只计算系统侧短路电流忽略逆变器型DG在故障瞬间的短路贡献可达额定电流1.2–2.0倍。解决在short_circuit_current函数中对每个DG节点叠加其短路电流function I_sc short_circuit_current(net, dg_site, cap_kW) % 系统侧短路电流简化为无穷大电源 I_sc_sys net.baseMVA / (sqrt(3)*net.baseKV*0.1); % 假设系统阻抗0.1pu % DG短路电流按1.5倍额定电流估算 I_dg_rated (cap_kW/1000) / (sqrt(3)*net.baseKV); % 额定电流(kA) I_sc_dg 1.5 * I_dg_rated; % 总短路电流 系统侧 DG侧近似相加 I_sc I_sc_sys I_sc_dg; end4.5 现象MATLAB运行报错“Out of memory”10节点网架都无法计算原因scan_capacity中未清理内存。每次调用潮流函数都会生成临时变量嵌套循环中不断累积尤其V,Iline等大矩阵未及时清除。解决在scan_capacity循环体内末尾添加% 循环内内存清理关键 clear V Iline Sloss Pdg Qdg; % 或更彻底使用parfor时用IdleTimeout控制同时在脚本开头设置% 启用内存优化 feature(Accelerator,on); % 减少图形缓存 set(0,DefaultFigureVisible,off);5. 进阶技巧用承载力灵敏度矩阵指导规划决策让评估结果直接驱动投资承载力评估的终极价值不是输出一份静态报告而是生成一张可操作的规划导航图。所谓“灵敏度”是指“单位容量DG接入对各约束指标的边际影响”。它能把抽象的“承载力低”转化为具体的“应优先改造哪条线路、加装何种设备”。我们通过有限差分法构建灵敏度矩阵并用它驱动三类关键决策。5.1 构建节点级灵敏度矩阵量化“接1kW光伏电压会降多少”对每个可接入节点i计算其对5个约束的灵敏度形成5×1向量s_ifunction sens_matrix compute_sensitivity(net, dg_site, dg_type, delta_cap) % delta_cap扰动量建议取5kW太小噪声大太大非线性失真 % 基准工况 [V0, ~] backward_forward_sweep(net, zeros(net.nbus,1), zeros(net.nbus,1)); base_Vmin min(abs(V0)); base_Vmax max(abs(V0)); base_Imax max_line_loading(net, V0, zeros(net.nbus,1), zeros(net.nbus,1)); base_Isc short_circuit_current(net, dg_site, 0); base_THD harmonic_distortion(net, dg_site, 0); % 扰动工况delta_cap Pdg zeros(net.nbus,1); Pdg(dg_site) delta_cap/1000; [V1, ~] backward_forward_sweep(net, Pdg, zeros(net.nbus,1)); pert_Vmin min(abs(V1)); pert_Vmax max(abs(V1)); pert_Imax max_line_loading(net, V1, Pdg, zeros(net.nbus,1)); pert_Isc short_circuit_current(net, dg_site, delta_cap); pert_THD harmonic_distortion(net, dg_site, delta_cap); % 灵敏度 (扰动后值 - 基准值) / delta_cap sens_matrix zeros(5,1); sens_matrix(1) (pert_Vmin - base_Vmin) / delta_cap; % 电压下限灵敏度 (pu/kW) sens_matrix(2) (pert_Vmax - base_Vmax) / delta_cap; % 电压上限灵敏度 (pu/kW) sens_matrix(3) (pert_Imax - base_Imax) / delta_cap; % 支路负载灵敏度 (pu/kW) sens_matrix(4) (pert_Isc - base_Isc) / delta_cap; % 短路电流灵敏度 (kA/kW) sens_matrix(5) (pert_THD - base_THD) / delta_cap; % 谐波灵敏度 (%/kW) end逻辑说明sens_matrix中每个元素代表“在该节点多接1kW DG对应约束指标的变化量”。例如sens_matrix(2)0.0015意味着接1kW光伏会使全网最高电压上升0.0015pu。若某节点sens_matrix(2)为正且绝对本文还有配套的精品资源点击获取