matlab编程潮流计算（请问，我用MATLAB利用PQ法做潮流计算时，为什么P和Q的误差开始越变越小，但迭代几次之后，会越变越大）

本文目录

• 请问，我用MATLAB利用PQ法做潮流计算时，为什么P和Q的误差开始越变越小，但迭代几次之后，会越变越大
• 如何用matlab分析潮流计算
• 请问用MATLAB编写潮流计算程序时，比如PQ法或者NR法，设置循环次数的意义在哪里呢
• MATLAB的电力系统潮流计算程序
• matlab用牛顿法计算潮流需要在命令窗口输入什么
• 用matlab进行牛顿拉夫逊法的潮流计算的编程，我也不会你能帮帮我吗万分感谢！！
• 求MATLAB简单潮流计算程序

请问，我用MATLAB利用PQ法做潮流计算时，为什么P和Q的误差开始越变越小，但迭代几次之后，会越变越大

1.首先要告诉你，不管是何种程序，在正确之前必须经过调试，因为错误是感觉不到的，必须经过程序逻辑的跟踪才能发现。当然，调试的过程一般比较辛苦，也需要一定的功力，建议可以请教一下比你水平高的人。否则，即便一个细小的错误也可能导致程序的错误，而你却永远无法发现！您所提到的现象，我相信没有人能够立刻回答为什么。 因此，下面提到的几点是给你调试的几点思路，帮助您自己找到问题，纠正错误！2.首先必须明白的是：潮流计算，牛拉法一般是：6次左右收敛，PQ分解法是12次左右收敛。首先要确定是你读取的数据正确，这个可以找一个正确的潮流程序进行验证，也可以将真值作为初始值在你的潮流程序中进行验证。你提到的迭代次数给人感觉好像是PQ分解法。3.将每次迭代后的V和相角都输出出来进行观察，你说迭代9次后，失配量反而变大了。到底是哪个节点的失配量变大了？是P还是Q？该节点的V和相角是否已经离谱了？ 再然后是确定，B, B,,两个矩阵是否是正确的。总共不就64个元素，一一进行检验。根据你描述的现象，我敢确定的是迭代9次后，你程序的V和 相角 已经发生了离谱的变化，希望你自己观察。4.检验程序是否有错误的逻辑，V和相角的修正量是+还是-，变压器支路的处理是否正确，等等。程序急不来，不能眼高手低。不经过深入的分析，不经过痛苦的修炼，永远不能成为高手，潮流计算也只是最最简单的入门。若是你的心态仅仅是为了应付一次毕业设计之类，一辈子再也不会做与程序有关的事情，我建议你找一位要好的老师或者同学来帮你去应付这个过程。

如何用matlab分析潮流计算

一、课程设计的目的与要求  目的：培养学生的电力系统潮流计算机编程能力，掌握计算机潮流计算的相关知识 要求：基本要求：  1.编写潮流计算程序； 2.在计算机上调试通过； 3.运行程序并计算出正确结果； 4.写出课程设计报告  二、设计步骤：  1.根据给定的参数或工程具体要求（如图），收集和查阅资料；学习相关软件（软件自选：本设计选择Matlab进行设计）。       2.在给定的电力网络上画出等值电路图。 3.运用计算机进行潮流计算。 4.编写设计说明书。 三、设计原理  1． 牛顿-拉夫逊原理  牛顿迭代法是取x0 之后，在这个基础上，找到比x0 更接近的方程的跟，一步一步迭代，从而找到更接近方程根的近似跟。牛顿迭代法是求方程根的重要方法之一，其最大优点是在方程f(x) = 0 的单根附近具有平方收敛，而且该法还可以用来求方程的重根、复根。电力系统潮流计算，一般来说，各个母线所供负荷的功率是已知的，各个节点电压是未知的（平衡节点外）可以根据网络结构形成节点导纳矩阵，然后由节点导纳矩阵列写功率方程，由于功率方程里功率是已知的，电压的幅值和相角是未知的，这样潮流计算的问题就转化为求解非线性方程组的问题了。为了便于用迭代法解方程组，需要将上述功率方程改写成功率平衡方程，并对功率平衡方程求偏导，得出对应的雅可比矩阵，给未知节点赋电压初值，一般为额课程设计报告  2  定电压，将初值带入功率平衡方程，得到功率不平衡量，这样由功率不平衡量、雅可比矩阵、节点电压不平衡量（未知的）构成了误差方程，解误差方程，得到节点电压不平衡量，节点电压加上节点电压不平衡量构成新的节点电压初值，将新的初值带入原来的功率平衡方程，并重新形成雅可比矩阵，然后计算新的电压不平衡量，这样不断迭代，不断修正，一般迭代三到五次就能收敛。  牛顿—拉夫逊迭代法的一般步骤： （1）形成各节点导纳矩阵Y。  （2）设个节点电压的初始值U和相角初始值e 还有迭代次数初值为0。 （3）计算各个节点的功率不平衡量。  （4）根据收敛条件判断是否满足，若不满足则向下进行。 （5）计算雅可比矩阵中的各元素。 （6）修正方程式个节点电压  （7）利用新值自第（3）步开始进入下一次迭代，直至达到精度退出循环。 （8）计算平衡节点输出功率和各线路功率 2． 网络节点的优化  1）静态地按最少出线支路数编号  这种方法由称为静态优化法。在编号以前。首先统计电力网络个节点的出线支路数，然后，按 出线支路数有少到多的节点顺序编号。当由n 个节点的出线支路相同时，则可以按任意次序对这n 个节点进行编号。这种编号方法的根据是导纳矩阵中，出线支路数最少的节点所对应的行中非零元素也  2）动态地按增加出线支路数最少编号在上述的方法中，各节点的出线支路数是按原始网络统计出 来的，在编号过程中认为固定不变的，事实上，在节点消去过程中，每消去一个节点以后，与该节点相连的各节点的出线支路数将发生变化(增加，减少或保持不变)。因此，如果每消去一个节点后，立即修正尚未编号节点的出线支路数，然后选其中支路数最少的一个节点进行编号，就可以预期得到更好的效果，动态按最少出线支路数编号方法的特点就是按出线最少原则编号时考虑了消去过程中各节点出线支路数目的变动情况。  3． MATLAB编程应用  Matlab 是“Matrix Laboratory”的缩写，主要包括：一般数值分析，矩阵运算、数字信号处理、建模、系统控制、优化和图形显示等应用程序。由于使用Matlab 编程运算与人进行科学计算的思路和表达方式完全一致，所以不像学习高级语言那样难于掌握，而且编程效率和计算效率极高，还可在计算机上直接输出结果和精美的图形拷贝，所以它的确为一高效的科研助手。

请问用MATLAB编写潮流计算程序时，比如PQ法或者NR法，设置循环次数的意义在哪里呢

在潮流计算程序中设这迭代次数上限是为了应对潮流发散的情况。比如一个正确的潮流计算 一般用牛拉法是 迭代6次，PQ分解法是 迭代12次。若潮流出现发散情况，迭代到20次，dP和dQ仍然不能收敛到你要的精度（比如0.0001），此时我们判定潮流发散，若不设一个迭代上限，可想而知，程序会一直迭代下去，此时数值已完全失去意义。 潮流发散意味着计算失败，这可能由两个原因导致：1数据存在问题，初值不在收敛域以内；2程序本身问题，考虑不全面，有bug。同时也可能产生两个现象，一种是 潮流直接发散，若干次迭代后，数值已经很离谱，此时数值肯定没有意义。还有一种是，雅克比矩阵接近病态，收敛极其缓慢（也可能是初值选取不当），此时适当的放大收敛次数限制，有一定可能得到收敛后的结果。多编程序，多思考，希望对你有帮助！

MATLAB的电力系统潮流计算程序

% 配电网前推回代潮流计算程序 % 使用IEEE 33节点配电系统作为算例，可实现弱环网情况下的潮流计算countnum=0;BranchData = ;                                      % 节点，负荷UB = 12.66;                             % 电压基准 kVSB = 10;                                % 功率基准 MVAZB = UB^2/SB;                           % 阻抗基准 ohmBranchData(:,) / ZB;     % 阻抗标幺化NodeData(:,) / SB / 1000;% 功率标幺化NN = 33;                                % 节点数A0 = zeros(NN);for n = 1:NN-1    A0(BranchData(n,1),BranchData(n,2)) = 1;end                                     % 形成 A0AssociatedMatrix=0;for n=2:NN-1    AssociatedMatrix(n,n)=1;    temp=BranchData(n-1,1);    AssociatedMatrix(n,1:n-1)=AssociatedMatrix(temp,1:n-1);endA0T = A0’;                              % 形成 A0 的转置S = ;       % 形成 SZL = ;       % 形成 ZLV = ones(NN,1);V(1) = 1;                            % 各个节点电压赋初值IL(NN,1) = -conj(S(NN) / V(NN));          % 最末支路电流赋初值Delta = 1;                              % 收敛判据赋初值TempV = V;                              % 赋初值，用于记忆上次迭代结果while Delta 》 1e-8    countnum=countnum+1;    IN = conj(S ./ V);                    % 节点注入电流        for n = 1:NN-1        IL(NN-n) = A0(NN-n,NN-n+1:end) * IL(NN-n+1:end) - IN(NN-n);    end                                 % 电流回代过程        for n = 2:NN        V(n) = A0T(n,1:n-1) * V(1:n-1) - ZL(n) * IL(n);    end                                 % 电压前推过程        Delta = max(abs(V-TempV));          % 更新收敛判据    TempV = V;                          % 记忆迭代结果endVangle(:,1)=abs(V);Vangle(:,2)=angle(V)/3.1415*180;for i=1:NN-1    st=BranchData(i,1);    en=BranchData(i,2);    Sij(i,1)=V(st)*conj((V(st)-V(en))/ZL(i+1));    Sji(i,1)=V(en)*conj((V(en)-V(st))/ZL(i+1));end

希望你能用到

matlab用牛顿法计算潮流需要在命令窗口输入什么

在公式（18）中，和分别表示状态变量与其修正量组成的列向量；为方阵，一般叫作雅可比矩阵，第i行j列元素为 ，它的大小为第i个函数对第j个变量求偏导；k则表示阵元素都在处取；同时，F(X)是由n个函数组成的n维列向量；在极坐标下，节点电压可如下表示：在这里插入图片描述（19）若和为已知大小的功率，与从节点电压求得的有功和无功功率之差，为功率的不平衡量，则节点功率不平衡量可用如下公式计算：在这里插入图片描述（20）节点功率可用各节点电压模值与相位表示，如下公式所示：在这里插入图片描述（21）式（21）中，为节点i和j的相位差。由以公式（18）-（21）推得牛顿法下，其潮流计算方程可写为：在这里插入图片描述（22）公式（22）中，雅可比矩阵的各元素为在这里插入图片描述（23）（24）（25）（26）在这里插入图片描述（27）（28）（29）（30）其中，节点导纳矩阵的元素由Gij 、Bij表示。随着国内外配电系统自动化水平不断提高，电力行业人员也开始更加深入地研究配电网系统。配电网潮流计算作为DMS（配电管理系统）的重要基础，受到广大行业界人士的关注。因此，配电网潮流计算，已然成为配电网分析的重要内容。配电网与输电网相比，两者有明显不同，前者一般采用网格结构，线路参数R/X的值较大，三相负荷不对称程度明显。这些特点使得在输电网中计算有效，如牛顿法，不再适用于配电网。为此，有学者提出了适用于配电网的潮流算法，主要包括基于回路方程的潮流算法、前推回推法和改进的牛顿-拉夫逊法（简称改进的牛拉法）。其中，基于回路方程的方法具有较强的网格处理能力和良好的收敛性，但该方法的节点数和分支数处理非常复杂。前推回推法是针对配电网的树状特性，可以避免潮流计算中的病态条件，同时速度更快。然而，由于其公式和算法与牛顿潮流算法不同，其在其它方面（如潮流优化）的应用将受到限制。改进牛顿法通过对传统法进行一定的近似，将J阵写成UDUT 的形式。U仅由网络拓扑决定，是一个上三角矩阵；D是一个对角矩阵。在牛拉法中，需要对J阵因子分解与前代回代，改进法则只有前推回代的计算过程。它很好地改善了传统法以及前推回推法。经过算例计算结果证明，改进法可以避免J阵病态，且拥有前推回代法的收敛速度、精度，又由于它属于牛顿型算法，所以该算法已经得到了广泛的运用。下面附带电力系统分析牛顿法算例及matlab程序：网络结构如下：系统结构图系统参数如下：在上图所示的简单电力系统中，系统中节点1、2为PQ节点，节点3为PV节点，节点4为平衡节点，已给定P1s+jQ1s=-0.30-j0.18 P2s+jQ2s=-0.55-j0.13 P3s=0.5 V3s=1.10 V4s=1.05∠0°容许误差ε=10-5节点导纳矩阵：导纳矩阵各节点电压：节点 e f v ζ1.0.984637 -0.008596 0.984675 -0.5001722.0.958690 -0.108387 0.964798 -6.4503063.1.092415 0.128955 1.100000 6.7323474.1.050000 0.000000 1.050000 0.000000各节点功率：节点 P Q1-0.300000 -0.1800002–0.550000 -0.1300003 0.500000 -0.5513054 0.367883 0.264698matlab程序如下：// 牛顿法潮流计算matlab程序clc;Y=;%导纳矩阵e=;%初始电压f=zeros(4,1);V=zeros(4,1);%节点电压Ws=;%初始功率W=zeros(6,1);n=length(Y);%节点数J=zeros(2*(n-1));%雅可比矩阵delta_v=zeros(1,6);delta_w=Ws;G=real(Y);B=imag(Y);S=zeros(4,2);c=0;%循环次数m=input(’请输入PQ节点数：’); while max(abs(delta_w))》10^-5for i=1:(n-1)%以下为求取雅可比矩阵 for j=1:(n-1) if (i~=j) J(2*i-1,2*j-1)=-(G(i,j)*e(i)+B(i,j)*f(i)); J(2*i,2*j)=-J(2*i-1,2*j-1); J(2*i-1,2*j)=B(i,j)*e(i)-G(i,j)*f(i); J(2*i,2*j-1)=J(2*i-1,2*j); end end endfor j=1:(n-2) J(6,2*j-1)=0; J(6,2*j)=0; end%以上为非对角线元素s1=0;s2=0;for i=1:(n-1) for j=1:n s1=s1+(G(i,j).*e(j)-B(i,j).*f(j)); s2=s2+(G(i,j).*f(j)+B(i,j).*e(j)); end J(2*i-1,2*i-1)=-s1-G(i,i) *e(i)-B(i,i)*f(i); J(2*i-1,2*i)=-s2+B(i,i) *e(i)-G(i,i)*f(i); s1=0; s2=0;endfor i=1:m for j=1:n s1=s1+G(i,j).*f(j)+B(i,j).*e(j); s2=s2+(G(i,j).*e(j)-B(i,j).*f(j)); end J(2*i,2*i-1)=s1+B(i,i) *e(i)-G(i,i)*f(i); J(2*i,2*i)=-s2+G(i,i) *e(i)+B(i,i)*f(i); s1=0; s2=0;endJ(6,5)=-2*e(3);J(6,6)=-2*f(3);%对角线元素求解for i=1:m for j=1:n s1=s1+e(i)*(G(i,j).*e(j)-B(i,j).*f(j))+f(i)*(G(i,j).*f(j)+B(i,j).*e(j)); s2=s2+f(i)*(G(i,j).*e(j)-B(i,j).*f(j))-e(i)*(G(i,j).*f(j)+B(i,j).*e(j)); end delta_w(2*i-1)=Ws(2*i-1)-s1; delta_w(2*i)=Ws(2*i)-s2; W(2*i-1)=s1; W(2*i)=s2; s1=0; s2=0;endfor j=1:n s1=s1+e(3)*(G(3,j).*e(j)-B(3,j).*f(j))+f(3)*(G(3,j).*f(j)+B(3,j).*e(j));enddelta_w(5)=Ws(5)-s1;delta_w(6)=(Ws(6)^2-(e(3)^2+f(3)^2));W(5)=s1;W(6)=sqrt(e(3)^2+f(3)^2);%以上求功率差值delta_v=-inv(J)*delta_w;for i=1:(n-1) e(i)=e(i)+delta_v(2*i-1); f(i)=f(i)+delta_v(2*i);end%求电压差值c=c+1; end for x=1:4 V(x)=e(x)+f(x)*1i; end%节点电压 s1=0; for x=3:4 for j=1:4 s1=s1+conj(Y(x,j))*conj(V(j)); end S(x,1)=real(V(x)*s1); S(x,2)=imag(V(x)*s1); s1=0; end%PV与平衡节点功率 for x=1:2 S(x,1)=W(2*x-1); S(x,2)=W(2*x); end%节点功率c JVS123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106运行结果如下：潮流计算结果文章知识点与官方知识档案匹配算法技能树首页概览34356 人正在系统学习中电磁MATLAB微信公众号程序从vx公众号获取，助力科研！打开CSDN APP，看更多技术内容matlab潮流程序,潮流计算的MATLAB源程序_weixin_39760433的博客-CSDN...1、简单潮流计算的牛顿拉夫逊程序,相关的原始数据数据数据输入格式如下:%B1是支路参数矩阵,第一列和第二列是节点编号。节点编号由小到大编写%对于含有变压器的支路,第一列为低压侧节点编号,第二列为高压侧节点%编号,将变压器的串联阻抗...继续访问...编程实现_某只旺仔的博客_电力系统潮流计算程序matlab潮流计算是电力系统中最基本,应用最广泛的一种计算,是电力系统稳定计算和故障分析的基础。这个代码通过matlab实现了用极坐标表示的牛顿法和P-Q分解法进行潮流计算,使用IEEE14节点系统进行测试,计算结果和应用matpower的潮流计算完全一致。代码...继续访问《em》matlab《/em》导数计算《p》matlab导数计算，通过编程实战掌握具体应用。包括matlab导数计算前、matlab导数计算中、matlab导数计算后。《/p》继续访问(完整版)基于MATLAB牛顿拉夫逊法进行潮流计算.doc》 %本程序的功能是用牛顿拉夫逊法进行潮流计算 n=input’ 请输入节点数 :n=; nl=input’ 请输入支路数 :nl=; isb=input’ 请输入平衡母线节点号 :isb=; pr=input’ 请输入误差精度 :pr=; B1=input’ 请输入由各支路参数形成的矩阵 :B1=; B2=input’ 请输入各节点参数形成的矩阵 :B2=; Y=zeros(n; e=zeros牛顿拉夫逊潮流计算matlab程序基于牛顿拉夫逊进行潮流计算，求得各节点电压，各支路功率流动，内附程序输入说明以及案例。珍藏多年的matlab潮流计算程序源代码集合，包含多个潮流计算程序【达摩老生出品，必属精品，亲测校正，质量保证】 资源名：珍藏多年的matlab潮流计算程序源代码集合，包含多个潮流计算程序 资源类型：matlab项目全套源码 源码说明： 全部项目源码都是经过测试校正后百分百成功运行的，如果您下载后不能运行可联系我进行指导或者更换。 适合人群：新手及有一定经验的开发人员牛顿-拉夫逊法潮流计算——采用极坐标形式采用极坐标形式的Newton-Raphson法进行潮流计算，并提供IEEE 14、57节点及New England 39节点的计算数据。使用方法：运行/修改PowerFlow_Newton_Polar.m文件。电力系统牛拉法潮流计算MATLAB程序实现潮流计算的MATLAB代码，使用牛顿拉夫逊法，可更改结点数，支路数和导纳阵，通用性强【课设/毕业设计】电力系统潮流计算（Matlab代码实现）潮流计算是电力系统分析中的一种最基本的计算，对给定系统进行潮流计算可以得到各母线上的电压、网络中的功率分布及功率损耗等。本文介绍了潮流计算在电力系统分析中的作用及潮流计算的发展状况。通过对潮流计算所用的数学模型进行分析，建立潮流计算的基本方程。牛顿-拉夫逊法是目前广泛应用的一种潮流计算方法，本文阐述了牛顿-拉夫逊潮流计算的基本原理。详细展示了利用matlab平台编写潮流计算程序的具体过程，通过实例证明基于牛顿-拉夫逊法的潮流计算程序具有收敛速度快、占用内存小的优点。潮流计算的程序实现手段。继续访问传统牛顿拉夫逊法潮流计算matlab程序.zip采用何仰赞《电力系统分析》中介绍的牛顿拉夫逊法matlab编程，程序还能实现N-1校核和线路网损分析，适合学习电气工程领域的学子。基于matpower的电力系统潮流计算matlab程序它是由美国康奈尔大学电力系统工程研究中心（PSERC of Cornell University）的RAY D. Zimmerman、Carlos E. Murillo-Sánchez和甘德强在Robert J. Thomas的指导下开发出来的，本章介绍的是MATPOWER4.0。每一个电网用变量名为“mpc”的结构体（structures）来定义，结构体mpc的不同字段用baseMVA、bus、branch、gen等来定义和返回电网的具体参数。列的数据类似于标准的IEEE 和PTI 列的数据格式。继续访问直流潮流计算matlab程序直流潮流发的特点是用电力系统的交流潮流（有功功率和无功功率）等值的直流电流来代替。甚至只用直流电路的解析法来分析电力系统的有功潮流，而不考虑无功分布对有功的影响。这样一来计算速度加快，但计算的准确度有所降低，本方法适用于对潮流计算准确度要求不高的计算场景。θ为网络中各节点的电压相位角的向量；P为节点注入的有功功率向量​。这就相当于线路两端的直流电位分别为θi和θj。（2）按照标幺值计算时，节点电压与其额定电压相差不大，故有：Ui≈Uj≈1.0；以IEEE9节点系统为算例，系统参数如下​。继续访问牛顿-拉夫逊法潮流计算matlab程序,牛顿—拉夫逊法潮流计算MATLAB程序.doc！！！！！！！！！！！！########################牛顿—拉夫逊法潮流计算程序By Yuluo%牛顿--拉夫逊法进行潮流计算n=input(’请输入节点数：n=’);n1=input(’请输入支路数：n1=’);isb=input(’请输入平衡母线节点号：isb=’);pr=input(’请输入误差精度：pr=’);B1=input(’请输入由支路参数形成的矩阵：B1=’);B...继续访问matlab程序 潮流计算,潮流计算matlab程序《潮流计算matlab程序》由会员分享，可在线阅读，更多相关《潮流计算matlab程序(3页珍藏版)》请在人人文库网上搜索。1、clear;%各节点参数：节点编号，类型，电压幅值，电压相位，注入有功，注入无功 %类型：1PQ节点，2PV节点，3平衡节点%本程序中将最后一个节点设为平衡节点R_1=1 1 1.0 0 0.2 0.2j;2 1 1.0 0 -0.45 -0.15j;3 1 1.0 0 ...继续访问matlab潮流计算函数,基于牛顿拉夫逊法潮流计算的matlab实验报告(含源程序和结果)...《基于牛顿拉夫逊法潮流计算的matlab实验报告(含源程序和结果)》由会员分享，可在线阅读，更多相关《基于牛顿拉夫逊法潮流计算的matlab实验报告(含源程序和结果)(24页珍藏版)》请在人人文库网上搜索。1、基于牛顿拉夫逊法潮流计算的matlab实验报告一、 实验目的和要求1.学习掌握matlab的基本用法2.应用MATLAB语言编写具有一定通用性的牛顿-拉夫逊法潮流计算程序。要求：(1)潮流计...继续访问matlab的做潮流计算,Matlab实现潮流计算程序关于Matlab趋势计算program.doc文件的实现，爱文共享的信息具有丰富的相关文档，每天都有成千上万的行业名人在该站点上共享最新信息.程序代码如下: 读取数据clcclearfilename =’txt’a = textread(filename)n = a(，)pinghengjd = a(，)phjddianya = a(，)jingdu = a(，)b = zeros (，)j = ...继续访问天然气潮流计算matlab程序此外，由于天然气在管道内传输时会产生压力损耗，所以每隔一段都要设置压缩机来提升节点压力，保证天然气系统的正常运行。（3）如果节点中存在压缩机，则可以求得压缩机消耗的等效天然气流量，将其作为一个负荷加入到潮流计算中；管道或者节点的约束，例如流量、节点压力、压缩机消耗电能和压缩比等，如下面四式所示。（2）根据天然气节点压力与节点间流量的关系式求出前一个节点的节点压力；（1）由用户负荷求得微型燃气轮机的功率????g，求出其消耗的天然气流量；（5）重复步骤（2）-（4），直至求解出整个天然气网络潮流。继续访问潮流计算程序matlab 牛拉法***隐藏网址***继续访问最新发布 基于matlab电力系统的牛顿-拉夫逊法进行潮流计算习题的具体计算继续访问热门推荐 潮流计算的matlab程序实现方法这是一个电气狗熬两个礼拜图书馆的成果，根据华中科技大学《电力系统分析》中原理编写，可用牛顿-拉夫逊和PQ分解法计算给定标幺值条件的潮流。本人水平有限，仅供参考，欢迎一起找Bug。 2019/11/17 添加算例系统图和基础数据、参考文献。 2019/01/05 添加word文档潮流计算课程设计。 2018/07/06 说明：由于本人变压器建模与PSASP不同，本人使用模型如下图，参数输入时...继续访问matlab潮流计算程序.rar_chaoliujisuan _matlab潮流计算_潮流程序_潮流计算_潮流计算程序在MATLAB环境下的潮流计算程序简单清晰【牛顿迭代】使用MATLAB综合设计电力系统潮流计算程序综合设计潮流计算程序，以数据文件输入电网数据原始数据，牛拉法算法实现，研究分析收敛特性及其与线路参数的关系。继续访问潮流计算matlab程序算法matlab写评论6130

用matlab进行牛顿拉夫逊法的潮流计算的编程，我也不会你能帮帮我吗万分感谢！！

B2=input(’请输入各节点参数形成的矩阵:B2=’);%本程序的功能是用牛顿-拉夫逊法进行潮流计算n=input(’请输入节点数:n=’);nl=input(’请输入支路数:nl=’);isb=input(’请输入平衡母线节电号:isb=’);pr=input(’请输入误差精度:pr=’);B1=input(’请输入由支路参数形成的矩阵:B1=’);%变压器侧为1，否则为0B2=input(’请输入各节点参数形成的矩阵:B2=’);X=input(’请输入由节点号及其对地阻抗形成的矩阵:X=’);X=input(’请输入由节点号及其对地阻抗形成的矩阵:X=’);Y=zeros(n);U=zeros(1,n);cta=zeros(1,n);V=zeros(1,n);O=zeros(1,n);S1=zeros(nl);for i=1:n if X(i,2)~=0; p=X(i,1); Y(p,p)=X(i,2); endendfor i=1:nl if B1(i,6)==0 p=B1(i,1);q=B1(i,2); else p=B1(i,2);q=B1(i,1); end Y(p,q)=Y(p,q)-1./(B1(i,3)*B1(i,5)); Y(q,p)=Y(p,q); Y(q,q)=Y(q,q)+1./(B1(i,3)*B1(i,5)^2)+B1(i,4)./2; Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4)./2;end %求导纳矩阵G=real(Y);B=imag(Y);for i=1:n cta(i)=angle(B2(i,3)); U(i)=abs(B2(i,3)); %V(i)=B2(i,4);endfor i=1:n S(i)=B2(i,1)-B2(i,2); B(i,i)=B(i,i)+B2(i,5);endP=real(S);Q=imag(S);ICT1=0;IT2=1;while IT2~=0 IT2=0;t1=1;t2=1; for i=1:n if i~=isb C(i)=0; D(i)=0; for j1=1:n C(i)=C(i)+U(i)*U(j1)*(G(i,j1)*cos(cta(i)-cta(j1))+B(i,j1)*sin(cta(i)-cta(j1))); D(i)=D(i)+U(i)*U(j1)*(G(i,j1)*sin(cta(i)-cta(j1))-B(i,j1)*cos(cta(i)-cta(j1))); end DP(t1)=P(i)-C(i); t1=t1+1; if B2(i,6)==2 DQ(t2)=Q(i)-D(i); t2=t2+1; end end end t1=t1-1;t2=t2-1; DPQ=; %求DP,DQ for i=1:t1+t2 if abs(DPQ(i))》pr IT2=IT2+1; end end H=zeros(t1,t1);N=zeros(t1,t2);K=zeros(t2,t1);L=zeros(t2,t2); for i=1:t1 for j1=1:t1 if j1~=isb&j1~=i H(i,j1)=0-U(i)*U(j1)*(G(i,j1)*sin(cta(i)-cta(j1))-B(i,j1)*cos(cta(i)-cta(j1))); elseif j1~=isb&j1==i H(i,j1)=U(i)^2*B(i,j1)+D(i); end end end for i=1:t1 for j1=1:t2 if j1~=isb&j1~=i N(i,j1)=0-U(i)*U(j1)*(G(i,j1)*cos(cta(i)-cta(j1))+B(i,j1)*sin(cta(i)-cta(j1))); elseif j1~=isb&j1==i N(i,j1)=0-U(i)^2*G(i,j1)-C(i); end end end for i=1:t2 for j1=1:t1 if j1~=isb&j1~=i K(i,j1)= U(i)*U(j1)*(G(i,j1)*cos(cta(i)-cta(j1))+B(i,j1)*sin(cta(i)-cta(j1))); elseif j1~=isb&j1==i K(i,j1)=U(i)^2*G(i,j1)-C(i); end end end for i=1:t2 for j1=1:t2 if j1~=isb&j1~=i L(i,j1)=0-U(i)*U(j1)*(G(i,j1)*sin(cta(i)-cta(j1))-B(i,j1)*cos(cta(i)-cta(j1))); elseif j1~=isb&j1==i L(i,j1)=U(i)^2*B(i,j1)-D(i); end end end J=;%求雅可比矩阵 modify=-J\DPQ; Dcta=modify(,:); t3=U(:,); DU=diag(t3,0)*modify(,:); t4=1; for i=1:t1 if B2(i,6)~=1 cta(1,i)=cta(1,i)+Dcta(t4,1); t4=t4+1; end end t5=1; for i=1:t2 if B2(i,6)==2 U(1,i)=U(1,i)+DU(t5,1); t5=t5+1; end end ICT1=ICT1+1;end %修正原值for i=1:n UU(i)=U(i)*cos(cta(i))+1i*U(i)*sin(cta(i));endfor p=1:n c(p)=0; for q=1:n c(p)=c(p)+conj(Y(p,q))*conj(UU(q)); end s(p)=UU(p)*c(p);enddisp(’--------------------------------------------------------------------------------’);disp(’各节点电压U为（节点从小到大排列）:’);disp(UU);disp(’--------------------------------------------------------------------------------’);disp(’各节点电压相角为（节点从小到大排列）:’);disp(180*angle(UU)/pi);disp(’--------------------------------------------------------------------------------’);disp(’按公式计算全部线路功率，结果如下:’);for i=1:nl if B1(i,6)==0 p=B1(i,1);q=B1(i,2); else p=B1(i,2);q=B1(i,1); end Si(p,q)=UU(p)*(conj(UU(p))*conj(B1(i,4)./2)+(conj(UU(p)*B1(i,5))-conj(UU(q)))*conj(1./(B1(i,3)*B1(i,5))));%各条支路首端功率Si f=; disp(f);endfor i=1:nl if B1(i,6)==0 p=B1(i,1);q=B1(i,2); else p=B1(i,2);q=B1(i,1); end Sj(q,p)=UU(q)*(conj(UU(q))*conj(B1(i,4)./2)+(conj(UU(q)./B1(i,5))-conj(UU(p)))*conj(1./(B1(i,3)*B1(i,5))));%各条支路末端功率Sj f=; disp(f);enddisp(’--------------------------------------------------------------------------------’);disp(’各条支路的功率损耗DS为(顺序同您输入B1时一样):’);for i=1:nl if B1(i,6)==0 p=B1(i,1);q=B1(i,2); else p=B1(i,2);q=B1(i,1); end DS(i)=Si(p,q)+Sj(q,p);%各条支路功率损耗DS disp(DS(i));endSp=0;for i=1:n Sp=Sp+UU(isb)*conj(Y(isb,i))*conj(UU(i));enddisp(’平衡节点的功率:’);disp(Sp);------------在网上找到的一个，希望能有帮助。

求MATLAB简单潮流计算程序

function lianxuchaoliuclear;clc;n=9;%节点数；nl=9;%支路数；isb=1;%平衡节点号；pr=0.00001;%误差精度；b1=;%依次是支路首端；末端，支路阻抗；对地电纳；支路变比；折算到哪一侧标志（高压侧为1；低压侧为0）；b2=;%依次是节点的发电机功率；负荷功率；节点电压初值；PV节点电压V给定值；节点无功补偿设备容量；节点分类标号（平衡1；PQ2；PV3）；Y=zeros(n);%求导纳阵；for i=1:nlif b1(i,6)==0 p=b1(i,1);q=b1(i,2); else p=b1(i,2);q=b1(i,1); end Y(p,q)=Y(p,q)-1./(b1(i,3)*b1(i,5)); Y(q,p)=Y(p,q); Y(q,q)=Y(q,q)+1./(b1(i,3)*b1(i,5)^2)+b1(i,4)./2; Y(p,p)=Y(p,p)+1./b1(i,3)+b1(i,4)./2; end%disp(’系统的导纳阵为：’);%disp(Y);G=real(Y);B=imag(Y);for i=1:ne(i)=real(b2(i,3)); f(i)=imag(b2(i,3)); v(i)=b2(i,4); endfor i=1:nS(i)=b2(i,1)-b2(i,2); B(i,i)=B(i,i)+b2(i,5); endP=real(S);Q=imag(S);w=zeros(2*n,1);Jac=zeros(2*n);t=0;while t==0for i=1:n if b2(i,6)~=isb C=0;D=0; for j=1:n C=C+G(i,j)*e(j)-B(i,j)*f(j); D=D+G(i,j)*f(j)+B(i,j)*e(j); end if b2(i,6)==2%P，Q节点； w(2*i)=P(i)-e(i)*C-f(i)*D; w(2*i-1)=Q(i)-f(i)*C+e(i)*D; else if b2(i,6)==3%P，V节点； w(2*i)=P(i)-e(i)*C-f(i)*D; w(2*i-1)=v(i)^2-(e(i)^2+f(i)^2); end end else w(2*i-1)=0; w(2*i)=0; end end%disp(w);w1=w(3:2*n);for i=1:nfor j=1:n if b2(i,6)~=isb if b2(i,6)==2%P，Q节点； if j~=i Jac(2*i,2*j-1)=-1*(G(i,j)*e(i)+B(i,j)*f(i)); Jac(2*i-1,2*j)=(G(i,j)*e(i)+B(i,j)*f(i)); Jac(2*i,2*j)=B(i,j)*e(i)-G(i,j)*f(i); Jac(2*i-1,2*j-1)=B(i,j)*e(i)-G(i,j)*f(i); else if j==i m=0;h=0; for r=1:n m=m+G(i,r)*e(r)-B(i,r)*f(r); h=h+G(i,r)*f(r)+B(i,r)*e(r); end Jac(2*i,2*j-1)=-1*m-G(i,i)*e(i)-B(i,i)*f(i); Jac(2*i-1,2*j)=-1*m+G(i,i)*e(i)+B(i,i)*f(i); Jac(2*i,2*j)=-1*h+B(i,i)*e(i)-G(i,i)*f(i); Jac(2*i-1,2*j-1)=h+B(i,i)*e(i)-G(i,i)*f(i); end end else if b2(i,6)==3%P，V节点； if j~=i Jac(2*i,2*j-1)=-1*(G(i,j)*e(i)+B(i,j)*f(i)); Jac(2*i-1,2*j)=0; Jac(2*i,2*j)=B(i,j)*e(i)-G(i,j)*f(i); Jac(2*i-1,2*j-1)=0; else if j==i m=0;h=0; for r=1:n m=m+G(i,r)*e(r)-B(i,r)*f(r); h=h+G(i,r)*f(r)+B(i,r)*e(r); end Jac(2*i,2*j-1)=-1*m-G(i,i)*e(i)-B(i,i)*f(i); Jac(2*i-1,2*j)=-2*f(i); Jac(2*i,2*j)=-1*h+B(i,i)*e(i)-G(i,i)*f(i); Jac(2*i-1,2*j-1)=-2*e(i); end end end end else Jac(2*i-1,2*j-1)=0; Jac(2*i,2*j)=0; Jac(2*i-1,2*j)=0; Jac(2*i,2*j-1)=0; end end end%disp(Jac);Jac1=Jac(3:2*n,3:2*n);for k=1:2*n-2m=0; for i=k+1:2*n-2 m=Jac1(i,k)./Jac1(k,k); w1(i)=w1(i)-m*w1(k); for j=k+1:2*n-2 Jac1(i,j)=Jac1(i,j)-m*Jac1(k,j); end end endE(2*n-2)=w1(2*n-2)./Jac1(2*n-2,2*n-2);for i=2*n-3:-1:1c=0; for k=i+1:2*n-2 c=c+Jac1(i,k)*E(k); E(i)=(w1(i)-c)./Jac1(i,i); end end%disp(E);for i=1:n-1e(i+1)=e(i+1)-E(2*i-1); f(i+1)=f(i+1)-E(2*i); end%disp(e);%disp(f);for i=1:n-1b(2*i-1)=abs(E(2*i-1)); b(2*i)=abs(E(2*i)); endKB=max(b);%disp(KB);if KB《prt=1; elset=0; endend%disp(e);%disp(f);for i=1:nfz(i)=sqrt(e(i)^2+f(i)^2); enddisp(fz);JJ1=zeros(2*n-1,2*n-1);Jac1=zeros(2*n-2,2*n-2);ttt=0;ccc=1;%迭代次数；L=0.05;T=0;%L是步长，T是参数；while ttt==0%确定扩展矩阵for i=1:n-1if b2(i+1,6)==3 JJ1(2*i,2*n-1)=1.5*1.05; JJ1(2*i-1,2*n-1)=0; else if b2(i+1,6)==2 if b2(i+1,2)==0 JJ1(2*i-1,2*n-1)=0; JJ1(2*i,2*n-1)=0; else JJ1(2*i,2*n-1)=-real(b2(i+1,2)); JJ1(2*i-1,2*n-1)=-imag(b2(i+1,2)); end end end endfor i=1:nfor j=1:n if b2(i,6)~=isb if b2(i,6)==2% P,Q节点； if j~=i Jac(2*i,2*j-1)=-1*(G(i,j)*e(i)+B(i,j)*f(i)); Jac(2*i-1,2*j)=(G(i,j)*e(i)+B(i,j)*f(i)); Jac(2*i,2*j)=B(i,j)*e(i)-G(i,j)*f(i); Jac(2*i-1,2*j-1)=B(i,j)*e(i)-G(i,j)*f(i); else if j==i m=0;h=0; for r=1:n m=m+G(i,r)*e(r)-B(i,r)*f(r); h=h+G(i,r)*f(r)+B(i,r)*e(r); end Jac(2*i,2*j-1)=-1*m-G(i,i)*e(i)-B(i,i)*f(i); Jac(2*i-1,2*j)=-1*m+G(i,i)*e(i)+B(i,i)*f(i); Jac(2*i,2*j)=-1*h+B(i,i)*e(i)-G(i,i)*f(i); Jac(2*i-1,2*j-1)=h+B(i,i)*e(i)-G(i,i)*f(i); end end else if b2(i,6)==3%P,V节点， if j~=i Jac(2*i,2*j-1)=-1*(G(i,j)*e(i)+B(i,j)*f(i)); Jac(2*i-1,2*j)=0; Jac(2*i,2*j)=B(i,j)*e(i)-G(i,j)*f(i); Jac(2*i-1,2*j-1)=0; else if j==i m=0;h=0; for r=1:n m=m+G(i,r)*e(r)-B(i,r)*f(r); h=h+G(i,r)*f(r)+B(i,r)*e(r); end Jac(2*i,2*j-1)=-1*m-G(i,i)*e(i)-B(i,i)*f(i); Jac(2*i-1,2*j)=-2*f(i); Jac(2*i,2*j)=-1*h+B(i,i)*e(i)-G(i,i)*f(i); Jac(2*i-1,2*j-1)=-2*e(i); end end end end else Jac(2*i-1,2*j-1)=0; Jac(2*i,2*j)=0; Jac(2*i-1,2*j)=0; Jac(2*i,2*j-1)=0; end end endJac1=Jac(3:2*n,3:2*n);for i=1:2*n-2for j=1:2*n-2 JJ1(i,j)=Jac1(i,j); end endif ccc==1%以负荷为连续参数；for j=1:2*n-2 JJ1(2*n-1,j)=0; end JJ1(2*n-1,2*n-1)=1; w2=zeros(2*n-1,1); for i=1:2*n-2 w2(i,1)=0; end w2(2*n-1,1)=1; end%disp(JJ1);if ccc》1%以切向量中分量最大（绝对值最大）的状态变量选定为连续参数；for i=1:2*n-2 Jd(i)=abs(d(i)); end for i=1:2*n-3 if Jd(i)》=Jd(i+1) zd=Jd(i);ek=i; else if Jd(i)《=Jd(i+1) zd=Jd(i+1);ek=i+1; end end end for j=1:2*n-1 JJ1(2*n-1,j)=0; end JJ1(2*n-1,ek)=1; w2=zeros(2*n-1,1); for i=1:2*n-1 w2(i,1)=0; end if d(ek)》0 w2(ek,1)=1; else if d(ek)《0 w2(ek,1)=-1; end end endfor k=1:2*n-1m=0; for i=k+1:2*n-1 m=JJ1(i,k)./JJ1(k,k); w2(i)=w2(i)-m*w2(k); for j=k+1:2*n-1 JJ1(i,j)=JJ1(i,j)-m*JJ1(k,j); end end endd(2*n-1)=w2(2*n-1)./JJ1(2*n-1,2*n-1);for i=2*n-2:-1:1c=0; for k=i+1:2*n-1 c=c+JJ1(i,k)*d(k); d(i)=(w2(i)-c)./JJ1(i,i); end end%disp(d);for i=1:n-1e(i+1)=e(i+1)+L*d(2*i-1); f(i+1)=f(i+1)+L*d(2*i); endT=T+L*d(2*n-1);disp(d(2*n-1));%disp(e);%disp(f);%disp(T);%对预估的近似解进行校正；tt=0;while tt==0for i=1:n if b2(i,6)~=isb C=0;D=0; for j=1:n C=C+G(i,j)*e(j)-B(i,j)*f(j); D=D+G(i,j)*f(j)+B(i,j)*e(j); end if b2(i,6)==2% P,Q节点； if b2(i,2)~=0 wj(2*i)=P(i)-T*real(b2(i,2))-e(i)*C-f(i)*D; wj(2*i-1)=Q(i)-T*imag(b2(i,2))-f(i)*C+e(i)*D; else wj(2*i)=P(i)-e(i)*C-f(i)*D; wj(2*i-1)=Q(i)-f(i)*C+e(i)*D; end else if b2(i,6)==3%P,V节点； wj(2*i)=P(i)+T*1.5*1.05-e(i)*C-f(i)*D; wj(2*i-1)=0; end end else wj(2*i)=0; wj(2*i-1)=0; end end wj1=wj(3:2*n); for i=1:n for j=1:n if b2(i,6)~=isb if b2(i,6)==2%P,Q节点； if j~=i JJac(2*i,2*j-1)=-1*(G(i,j)*e(i)+B(i,j)*f(i)); JJac(2*i-1,2*j)=(G(i,j)*e(i)+B(i,j)*f(i)); JJac(2*i,2*j)=B(i,j)*e(i)-G(i,j)*f(i); JJac(2*i-1,2*j-1)=B(i,j)*e(i)-G(i,j)*f(i); else if j==i m=0;h=0; for r=1:n m=m+G(i,r)*e(r)-B(i,r)*f(r); h=h+G(i,r)*f(r)+B(i,r)*e(r); end JJac(2*i,2*j-1)=-1*m-G(i,i)*e(i)-B(i,i)*f(i); JJac(2*i-1,2*j)=-1*m+G(i,i)*e(i)+B(i,i)*f(i); JJac(2*i,2*j)=-1*h+B(i,i)*e(i)-G(i,i)*f(i); JJac(2*i-1,2*j-1)=h+B(i,i)*e(i)-G(i,i)*f(i); end end else if b2(i,6)==3%P,V节点； if j~=i JJac(2*i,2*j-1)=-1*(G(i,j)*e(i)+B(i,j)*f(i)); JJac(2*i-1,2*j)=0; JJac(2*i,2*j)=B(i,j)*e(i)-G(i,j)*f(i); JJac(2*i-1,2*j-1)=0; else if j==i m=0;h=0; for r=1:n m=m+G(i,r)*e(r)-B(i,r)*f(r); h=h+G(i,r)*f(r)+B(i,r)*e(r); end JJac(2*i,2*j-1)=-1*m-G(i,i)*e(i)-B(i,i)*f(i); JJac(2*i-1,2*j)=-2*f(i); JJac(2*i,2*j)=-1*h+B(i,i)*e(i)-G(i,i)*f(i); JJac(2*i-1,2*j-1)=-2*e(i); end end end end else JJac(2*i-1,2*j-1)=0; JJac(2*i,2*j)=0; JJac(2*i-1,2*j)=0; JJac(2*i,2*j-1)=0; end end end JJac1=JJac(3:2*n,3:2*n); for k=1:2*n-2 m=0; for i=k+1:2*n-2 m=JJac1(i,k)./JJac1(k,k); wj1(i)=wj1(i)-m*wj1(k); for j=k+1:2*n-2 JJac1(i,j)=JJac1(i,j)-m*JJac1(k,j); end end end E1(2*n-2)=wj1(2*n-2)./JJac1(2*n-2,2*n-2); for i=2*n-3:-1:1 c=0; for k=i+1:2*n-2 c=c+JJac1(i,k)*E1(k); E1(i)=(wj1(i)-c)./JJac1(i,i); end end %disp(E1); for i=1:n-1 e(i+1)=e(i+1)-E1(2*i-1); f(i+1)=f(i+1)-E1(2*i); end for i=1:n-1 bx(2*i-1)=abs(E1(2*i-1)); bx(2*i)=abs(E1(2*i)); end KB1=max(bx); if KB1《pr tt=1; else tt=0; end end%disp(e);%disp(f);if d(2*n-1)》0ttt=0;ccc=ccc+1; elsettt=1; endenddisp(T);disp(ccc);%disp(e);%disp(f);for i=1:nfz1(i)=sqrt(e(i)^2+f(i)^2); enddisp(fz1);