电力系统计算机辅助分析试验报告

1、附录：学生实验报告表头格式电力工程学院学生实验报告实验课程名称：电力系统计算机辅助分析开课实验室：计算中心4082014年 12月 11日年级、专业、班学号姓名成绩实验项目名电力系统计算机辅助分析指导称教师教师 评语教师签名：年月日实验一MATLAB 软件的基本操作及程序的基本结构一、程序流程、源码及仿真结果启动 MATLAB 软件，识别出常用的四个窗口命令窗口（Command Window ）、历史命令窗口（ Command History ）、工作空间浏览器（ Workspace Browser ）。熟悉 Matlab 工作环境。2、上机求解下面的电路问题，脚本文件已经给出。已知某个电路的

2、端电压的表达式为：u 6e 2 t , t0 10 区间时，试绘制电压的波形。提示：在MATLAB 的 M 文件编辑器中键入以下命令语句，并保存为 ex.m： Clear,clc,closet=0:0.5:10;u=6*exp(-2*t);plot(t,u );title(电路RC的电压响应曲线 );xlabel(时间 /s );ylabel(电压 /v)grid%RC电路的电压响应曲线654v/压 3电210012345678910时 间 /s5、 对于如下电路图，已知us (t)200cos( t 90 o ) V， R122， L50mH ， R233 ，C100 F ，自己编写程序，画

3、出 uc (t) 和 us (t)的波形（四个周期）。R1L+u s(t )Cu c ( t )R2-提示：写出 us(t ) 的向量表达式用于uc (t ) 的计算（写出串联段和并联段的阻抗，利用分压公式），写出us (t) 和 uc (t ) 的幅值和相角就能写出波形表达式，设定波形的长度（四个周波）和采样频率（两点间隔）后作图。注意画图时使用标题，横纵坐标及网格线等命令。3.%clear,clc,closemUsm=200;% 的幅值aUsm=(pi/180)*(+90);% 的相角Usm=mUsm*exp(j*aUsm);% 用于计算w=100*pi;r1=22;l=0.05;r2=3

4、3;c=0.0001;z1=r1+j*w*l;z2=r2/(j*w*c*r2+1);Ucm=Usm*z2/(z1+z2);mUcm=abs(Ucm)%求幅值aUcm=angle(Ucm)%求相角T=4*2*pi/w;%用 w 表示 ，T 表示四个周波N=100;%画图的点数dt=T/N;%两点间隔，即采样频率t=0:dt:T;for k=1:101%画波形us(k)=mUsm*cos(w*t(k)+aUsm);uc(k)=mUcm*cos(w*t(k)+aUcm);endplot(t,us,t,uc)%注意区分两条曲线的线型和颜色grid ontitle()xlabel(,fontsize,1

5、6)% fontsize ,16 表示横坐标的字体大小为 16ylabel()legend(us,uc)%一、 讨论1、 MATLAB软件是什么？有什么特点？答：Matlab 是一个高级的矩阵/ 阵列语言，它包含控制语句、函数、数据结构、输入和输出和面向对象编程特点。用户可以在命令窗口中将输入语句与执行命令同步，也可以先编写好一个较大的复杂的应用程（M 文件）后再一起运行。特点： 1. 编程效率高， 2.用户使用方便， 3. 扩充能力强，交互性好，4. 移植性和开放性很好， 5. 语句简单，内涵丰富，6. 高效方便的矩阵和数组运算， 7.方便的绘图功能。3.如何获得命令或函数的帮助文档？答：M

6、ATLAB 的各个函数，不管是内建函数、 M 文件函数、还是 MEX 文件函数等，一般它们都有说明，各个工具箱通常情况下也具有一个与工具箱名相间的 M 文件用来说明工具箱的构成内容等，M 文件的使用帮助和函数功能在 MATLAB 命令窗口中， 可以通过指令来获取这些纯文本的帮助信息。通常能够起到帮助作用，获取帮助信息的指令有help 、lookfor 、 which 、 doc、get 、type等。实验二电力系统计算中常用的数值算法及电力网络的数学模型1、用因子表法求解电路方程中的节点电压；A= 1/10+1/15 -1/10 -1/15 -1-1/10 1/10+1/50+1/30 -1/

7、50 0-1/15 -1/50 1/15+1/50+1/30 01000;B=0 0 0 10;n,m=size(A);for i=1:nA(i,i)=1/A(i,i);for j=i+1:nA(i,j)=A(i,j)*A(i,i);endfor k=i+1:nfor j=i+1:nA(k,j)=A(k,j)-A(k,i)*A(i,j);endendenddisp(矩阵 A 的因子表为： );disp(A)for i=1:nB(i)=B(i)*A(i,i);for j=i+1:nB(j)=B(j)-A(j,i)*B(i);endendfor i=n-1:-1:1for j=i+1:-1:2B(

8、j-1)=B(j-1)-A(j-1,i+1)*B(i+1);endenddisp(在因子表的基础上求解线性方程组的解为：x=);disp(B)2、将图 2.1 中的电压源分别改为2V、 4V、 6V、 8V 时，用因子表进行求解，体会该算法带来的益处；电压源改为2V电压源改为4V电压源改为6V电压源改为8V3、 上机验证课本例17.617.7 的计算结果；例 17.6例 17.74、 设计系统变更时修改导纳矩阵的函数（参考课本码。P.323），并调试通过，给出导纳矩阵修正函数的调用规则及其主要的实现代% n=input(请输入节点数 :n=);% nl=input(请输入支路数 :nl=);%

9、 B=input(请输入由支路参数形成的矩阵 :B=);% X=input(请输入由节点号及其对地阻抗形成的矩阵 :X=); n=4;nl=4;B=1 2 0.08+0.4i 0 1 0;2 3 0.1+0.4i 0 1 0;3 4 0.3i 0 1 0;1 3 0.12+0.5i 0 1 0; X=1 0;2 0;3 0;4 0Y=zeros(n);for i=1:nif X(i,2)=0;p=X(i,1);Y(p,p)=1./X(i,2);endendfor i=1:nlif B(i,6)=0p=B(i,1); q=B(i,2);elsep=B(i,2); q=B(i,1);endY(p,

10、q)=Y(p,q)-1./(B(i,3)*B(i,5);Y(q,p)=Y(p,q);Y(q,q)=Y(q,q)+1./(B(i,3)*B(i,5)2)+B(i,4)./2;Y(p,p)=Y(p,p)+1./B(i,3)+B(i,4)./2;enddisp(导纳矩阵 Y=:);disp(Y)B1=input(请输入追加支路参数形成的矩阵:B1=);% B1=3 4 0.3i 0 1.2 0;% B1=1 2 0.3i 0 1 0;% B1=1 4 0.3i 0 1 0;% B1=3 5 0.3i 0 1 0; z=B1(3);% if B1(6)=0% p=B1(1); q=B1(2);% el

11、se% p=B1(2); q=B1(1);% endp=B1(1); q=B1(2);if B1(6)=1B1(5)=1/(B1(5);endif B1(2)nY(p,p)=Y(p,p)+1/z; Y(p,q)=-1/z; Y(q,p)=-1/z; Y(q,q)=(1/z)+B1(4);endif Y(p,q)=0Y(p,p)=Y(p,p)+1/z;Y(q,q)=Y(q,q)+1/z;Y(p,q)=-1/z; Y(q,p)=-1/z;endif Y(p,q)=0&B1(2)=nfor i=1:nlif B(i,1)=B1(1)&B(i,2)=B1(2)if B(i,5)=B1(5)Y(p,p)

12、=Y(p,p)+1/z;Y(q,q)=Y(q,q)+1/z;Y(p,q)=Y(p,q)-1/z;Y(q,p)=Y(q,p)-1/z;elseK1=B1(5);K=B(i,5);Y(p,p)=Y(p,p)+(K12-K2)*(1/z);Y(q,q)=Y(q,q);Y(p,q)=Y(p,q)-(K1-K)*(1/z);Y(q,p)=Y(q,p)-(K1-K)*(1/z);endendendendifB1(2)nn=n+1;Y1=zeros(n);for i=1:nfor j=1:nY1(i,j)=Y(i,j);endendenddisp(新的导纳矩阵Y1=);disp(Y1);仿真结果：X =10

13、203040导纳矩阵 Y=:0.9346 - 4.2949i-0.4808 + 2.4038i-0.4539 + 1.8911i0-0.4808 + 2.4038i1.0690 - 4.7568i-0.5882 + 2.3529i0-0.4539 + 1.8911i-0.5882 + 2.3529i1.0421 - 7.5773i0 + 3.3333i000 + 3.3333i0 - 3.3333i请输入追加支路参数形成的矩阵 :B1=3 4 0.3i 0 1.2 0;新的导纳矩阵 Y1=0.9346 - 4.2949i-0.4808 + 2.4038i-0.4539 + 1.8911i0-0

14、.4808 + 2.4038i1.0690 - 4.7568i-0.5882 + 2.3529i0-0.4539 + 1.8911i-0.5882 + 2.3529i1.0421 - 9.0440i0 + 4.0000i000 + 4.0000i0 - 3.3333i思考题：1、高斯消去法与因子表法有何异同？答： 因子表法是以高斯消去法为基础求解线性方程组的一种方法，当线性方程组的系数矩阵不变，而常数项发生变化时，采用因子表法可避免系数矩阵的重复运算。因子表法首先要建立因子表，然后用因子表对常数项进行消去、回代，从而求出线性方程组的解。2、节点导纳阵与节点阻抗阵之间有何关系，二者之间如何转换？

15、电力系统的网络方程式一般都用节点方程式表示，节点方程有导纳型和阻抗型两种。节点导纳矩阵具有对称、稀疏、可以根据网络接线图直接写出、易于修改的特点；节点阻抗矩阵是节点导纳矩阵的逆阵，它具有对称、满阵的特点，它不能由网络接线图直接写出，通常由支路追加法形成。实验三短路电流计算合并课本第 18章例程。即要求实现计算指定节点各种故障条件下，各支路的故障电流分布情况。其中对于实用算法，可只提供计算电抗。提示：尽量使用先前实验中已封装好的函数。2.用新程序验证第18 章的各例题。（注意短路电流曲线绘制时，至少包含57 个点）W=3;switch Wcase1clcclear allNF=1for i1=1

16、:NFclear;n1=3n2=3n0=3nl1=5nl2=5nl0=5Lf=1%单相接地短路为1 横向故障为0f=3If=0zf=0zg=0B1=0 1 0.15i 0 1 0;0 2 0.075i 0 1 0;1 2 0.1i 0 1 0;1 3 0.1i 0 1 0;2 3 0.1i 0 1 0B2=0 1 0.15i 0 1 0;0 2 0.075i 0 1 0;1 2 0.1i 0 1 0;1 3 0.1i 0 1 0;2 3 0.1i 0 1 0B0=0 1 0.05i 0 1 0;0 2 0.025i 0 1 0;1 2 0.2i 0 1 0;1 3 0.2i 0 1 0;2 3

17、 0.2i 0 1 0A1=1 1 1;1 -1./2-i*sqrt(3)/ 2 -1./2+i*sqrt(3)/ 2;1 -1./2+i*sqrt(3)/ 2 -1./2-i*sqrt(3)/ 2Vcs=1;1;1Z1=zeros(n1);Z2=zeros(n2);Z0=zeros(n0);Y1=zeros(n1);Y2=zeros(n2);Y0=zeros(n0);V1=zeros(n1,1);V2=zeros(n2,1);V0=zeros(n0,1);I1=zeros(nl1,1);I2=zeros(nl2,1);I0=zeros(nl0,1);for m1=1:3m=0;if m1=1

18、nl=nl1;B=B1;Z=Z1;n=n1;elseif m1=2nl=nl2;B=B2;Z=Z2;n=n2;elseif m1=3nl=nl0;B=B0;Z=Z0;n=n0;endm=0;Z=zeros(n);V=zeros(n);I=zeros(nl);for k1=1:nlp=B(k1,1);q=B(k1,2);if B(k1,6)=0k=1./B(k1,5);else k=B(k1,5);endif p=0if qm%追加接地树支Z(q,q)=B(k1,3);m=m+1;else%追加接地连支for i=1:mZ(i,m+1)=-Z(i,q);Z(m+1,i)=-Z(q,i);endZ

19、(m+1,m+1)=Z(q,q)+B(k1,3);for i=1:mfor j=1:mZ(i,j)=Z(i,j)-Z(i,m+1)*Z(m+1,j)./Z(m+1,m+1);endZ(i,m+1)=0endfor i=1:m+1Z(m+1,i)=0endendelse if qm%追加不接地树支for i=1:mZ(i,q)=Z(i,p)*k;Z(q,i)=Z(p,i)*k;endZ(q,q)=k2*Z(p,p)+k2*B(k1,3);m=m+1;elsefor i=1:m%追加不接地连支Z(i,m+1)=k*Z(i,p)-Z(i,q);Z(m+1,i)=k*Z(p,i)-Z(q,i);end

20、Z(m+1,m+1)=k2*Z(p,p)+Z(q,q)-2*k*Z(p,q)+k2*B(k1,3);for i=1:mfor j=1:mZ(i,j)=Z(i,j)-Z(i,m+1)*Z(m+1,j)./Z(m+1,m+1);endZ(i,m+1)=0;endfor i=1:m+1Z(m+1,i)=0;endendendendif m1=1Z1=Z;elseif m1=2Z2=Z;elseif m1=3Z0=Z;endend%以上程序求出正 ,负,零序对应的阻抗阵if If=0%求出横向故障的阻抗阵Z1(f,f)=Z1(f,f);Z2(f,f)=Z2(f,f);Z0(f,f)=Z0(f,f);e

21、lse%求出纵向故障的阻抗阵Z1(f,f)=Z1(f,f)+Z1(1,1)-2*Z1(f,1);Z2(f,f)=Z2(f,f)+Z2(1,1)-2*Z2(f,1);Z0(f,f)=Z0(f,f)+Z0(1,1)-2*Z0(f,1);end%-故障点自阻抗 -ifLf=1%求出正负零序短路节点的电流标幺值Z4=Z2(f,f)+Z0(f,f)+3*zf;k2=1;k0=1;I1(f)=Vcs(f)./(Z1(f,f)+Z4);I2(f)=k2*I1(f);I0(f)=k0*I1(f);a1=I1(f);a2=I2(f);a0=I0(f);elseif Lf=2%两相短路a2=I2(f);a0=I0

22、(f);elseif Lf=3%两相接地短路Z4=zf+(Z2(f,f)+zf)*(Z0(f,f)+zf+3*zg)./(Z2(f,f)+Z0(f,f)+2*zf+3*zg);k2=-(Z0(f,f)+zf+3*zg)./(Z2(f,f)+Z0(f,f)+2*zf+3*zg); k0=-(Z2(f,f)+zf)./(Z2(f,f)+Z0(f,f)+2*zf+3*zg);I1(f)=Vcs(f)./(Z1(f,f)+Z4);I2(f)=k2*I1(f);I0(f)=k0*I1(f);a1=I1(f); a2=I2(f);a0=I0(f);elseif Lf=4%单相断线Z4=zf+(Z2(f,f

23、)+zf)*(Z0(f,f)+zf)./(Z2(f,f)+Z0(f,f)+2*zf);k2=-(Z0(f,f)+zf)./(Z2(f,f)+Z0(f,f)+2*zf); k0=-(Z2(f,f)+zf)./(Z2(f,f)+Z0(f,f)+2*zf);I1(f)=Vcs(f)./(Z1(f,f)+Z4);I2(f)=k2*I1(f);I0(f)=k0*I1(f);a1=I1(f);a2=I2(f);a0=I0(f);elseif Lf=5%两相断线Z4=Z2(f,f)+Z0(f,f)+3*zf;k2=1;k0=1;I1(f)=Vcs(f)./(Z1(f,f)+Z4);I2(f)=k2*I1(f

24、);I0(f)=k0*I1(f);a1=I1(f);a2=I2(f);a0=I0(f); enddisp(故障点处的正序电流标幺值为:);disp(a1);disp(故障点处的负序电流标幺值为:);disp(a2);disp(故障点处的零序电流标幺值为:);disp(a0);for i=1:n1%求出各节点的电压标幺值V1(i)=Vcs(i)-I1(f)*Z1(i,f);%正序endfor i=1:n2V2(i)=-I2(f)*Z2(i,f);%负序endfor i=1:n0V0(i)=-I0(f)*Z0(i,f);%零序enddisp(各节点的正序电压标幺值(节点号从小到大排):);disp

25、(V1);disp(各节点的负序电压标幺值(节点号从小到大排):);disp(V2);disp(各节点的零序电压标幺值(节点号从小到大排):);disp(V0);for i=1:3%求各支路的各序电流标幺值if i=1B=B1;V=V1;nl=nl1;elseif i=2B=B2;V=V2;nl=nl2;elseif i=3B=B0;V=V0;nl=nl0;endfor j=1:nlif B(j,6)=0k=B(j,5);else k=1./B(j,5);endp=B(j,1);q=B(j,2);if p=0e=0;b=B(j,3);I(j)=(e-V(q)./k)./b;%接地支路电流els

26、eI(j)=(V(p)-V(q)./k)./B(j,3);%不接地支路电流endendif i=1I1=I;elseif i=2I2=I;elseif i=3I0=I;endenddisp(各支路的正序电流标幺值为(顺序同您输入 B 时一样 ):);for i=1:nldisp(I1(i);enddisp(各支路的负序电流标幺值为(顺序同您输入 B 时一样 ):);for i=1:nldisp(I2(i);enddisp(各支路的零序电流标幺值为(顺序同您输入 B 时一样 ):);for i=1:nldisp(I0(i);endfor i=1:3%求故障点 A,B,C三相的总电流 ,各个节点

27、A,B,C三相各相的总电压 ,各支路 A,B,C三相各相的总电流if i=1%故障点电流s0=a0;s1=a1;s2=a2;elseif i=2%各序支路电压s0=V0;s1=V1;s2=V2;else%各序支路电流s0=I0;s1=I1;s2=I2;endfor j=1:3%不对称系数CA=s0*A1(j,1)+s1*A1(j,2)+s2*A1(j,3);if j=1A=CA;elseif j=2B=CA;elseif j=3C=CA;endendif i=1Iadl=A;Ibdl=B;Icdl=C;elseif i=2Vajy=A;Vbjy=B;Vcjy=C;elseIazl=A;Ibzl

28、=B;Iczl=C;endenddisp(故障点处 A 相电流 Iadl 为:);disp(Iadl);disp( 故障点处 B 相电流 Ibdl 为 :);disp(Ibdl);disp(故障点处 C 相电流 Icdl 为:);disp(Icdl);disp( 各节点处 A 相电压 Vajy 分别为 (节点号从小到大排 ):);disp(Vajy); disp(各节点处 B 相电压 Vbjy 分别为 (节点号从小到大排 ):)disp(Vbjy); disp(各节点处 C 相电压 Vcjy 分别为 (节点号从小到大排 ):);disp(Vcjy);disp(各支路 A 相电流 Iazl 分别

29、为 (顺序同你输入B 时一样 ):);for i=1:nldisp(Iazl(i);enddisp(各支路 B 相电流 Ibzl 分别为 (顺序同你输入B 时一样 ):);for i=1:nldisp(Ibzl(i);enddisp(各支路 C 相电流 Iczl 分别为 (顺序同你输入B 时一样 ):);for i=1:nldisp(Iczl(i);endendcase 2NF=1;n=4;nl=5;B=0 1 0.2i 0 1 0;0 2 4i 0 1 0;1 3 0.51i 0 1 0;2 3 0.59i 0 1 0;3 4 1.43i 0 1 0 ;D=4,0;V0=1;1;1;1; m

30、=0;Z=zeros(n);V=zeros(n);I=zeros(nl);for k1=1:nlp=B(k1,1);q=B(k1,2);if B(k1,6)=0k=1./B(k1,5);else k=B(k1,5);endif p=0if qmZ(q,q)=B(k1,3);m=m+1;elsefor i=1:m,Z(i,m+1)=-Z(i,q);Z(m+1,i)=-Z(q,i);endZ(m+1,m+1)=Z(q,q)+B(k1,3);for i=1:mfor j=1:mZ(i,j)=Z(i,j)-Z(i,m+1)*Z(m+1,j)./Z(m+1,m+1);endZ(i,m+1)=0;endf

31、or i=1:m+1Z(m+1,i)=0;endendelse if qmfor i=1:mZ(i,q)=Z(i,p)*k;Z(q,i)=Z(p,i)*k;endZ(q,q)=k2*Z(p,p)+k2*B(k1,3);m=m+1;elsefor i=1:mZ(i,m+1)=k*Z(i,p)-Z(i,q);Z(m+1,i)=k*Z(p,i)-Z(q,i);endZ(m+1,m+1)=k2*Z(p,p)+Z(q,q)-2*k*Z(p,q)+k2*B(k1,3);for i=1:mfor j=1:mZ(i,j)=Z(i,j)-Z(i,m+1)*Z(m+1,j)./Z(m+1,m+1);endZ(i,

32、m+1)=0;endfor i=1:m+1Z(m+1,i)=0;endendendendfor k=1:NFI(D(k,1),D(k,1)=V0(D(k,1),1)./(Z(D(k,1),D(k,1)+D(k,2);ft=num2str(D(k,1);ts1=(点短路时 );ts2=( 电流的标么值If=);dn=strcat(ft,ts1,ts2);disp(dn);disp(I(D(k,1),D(k,1);for i=1:nV(i,i)=V0(i,1)-I(D(k,1),D(k,1)*Z(i,D(k,1);endfor i=1:nlif B(i,6)=0K=B(i,5);else k=1.

33、/B(i,5);endp=B(i,1);q=B(i,2);if p=0e=0;b=B(i,3);I(i,i)=(e-V(q,q)./k)./b;elseI(i,i)=(e-V(p,p)-V(q,q)./k)./B(i,3);endenddisp(各节点的电压标么值U 为：);for i=1:ndisp(V(i,i);enddisp(各支路短路电流的标么值I 为： );for i=1:nldisp(I(i,i);endendcase 3clcNF=1;n=3;nl=5;Sb=100;D=3 0;vb=115 115 115;V0=1;1;1;B=0 1 0.15i 0 1 0 100 1 1;0

34、 2 0.075i 0 1 0 200 1 1;1 2 0.1i 0 1 0 0 0 2;1 3 0.1i 0 1 0 0 0 2;2 3 0.1i 0 1 0 0 0 2; N=5;m=0;Z=zeros(n);sq3=sqrt(3);s=zeros(N,1);for k1=1:nlp=B(k1,1);q=B(k1,2);if B(k1,6)=0k=1./B(k1,5);else k=B(k1,5);endif p=0if qmZ(q,q)=B(k1,3);m=m+1;elsefor i1=1:m,Z(i1,m+1)=-Z(i1,q);Z(m+1,i1)=-Z(q,i1);endZ(m+1,

35、m+1)=Z(q,q)+B(k1,3);for i1=1:mfor j=1:mZ(i1,j)=Z(i1,j)-Z(i1,m+1)*Z(m+1,j)./Z(m+1,m+1);endZ(i1,m+1)=0;endfor i1=1:m+1Z(m+1,i1)=0;endendelse if qmfor i1=1:mZ(i1,q)=Z(i1,p)*k;Z(q,i1)=Z(p,i1)*k;endZ(q,q)=k2*Z(p,p)+k2*B(k1,3);m=m+1;elsefor i1=1:mZ(i1,m+1)=k*Z(i1,p)-Z(i1,q);Z(m+1,i1)=k*Z(p,i1)-Z(q,i1);end

36、Z(m+1,m+1)=k2*Z(p,p)+Z(q,q)-2*k*Z(p,q)+k*B(k1,3);for i1=1:mfor j=1:mZ(i1,j)=Z(i1,j)-Z(i1,m+1)*Z(m+1,j)./Z(m+1,m+1);endZ(i1,m+1)=0;endfor i1=1:m+1Z(m+1,i1)=0;endendendendVb=zeros(1,n);V=zeros(1,n);Vd=zeros(1,n);Ib=zeros(1,nl);for r=1:NFIdb=V0(D(r,1),1)./(Z(D(r,1),D(r,1)+D(r,2);ts1=(以下是短路点 );ft=NUm2st

37、r(D(r,1);ts2=(各时刻的结果 );dn=strcat(ts1,ft,ts2);disp(dn);for k=1:nVb(k)=V0(k,1)-Idb*Z(k,D(r,1);endfor i1=1:nlp=B(i1,1);q=B(i1,2);if p=0&B(i1,8)=0if B(i1,6)=0k=B(i1,5);V(i1)=Vb(p)-Vb(q)./kelsek=1./B(i1,5);V(i1)=k*Vb(p)-Vb(q);endelseV(i1)=1-Vb(q);endIb(i1)=V(i1)./B(i1,3);enddisp(0s 时短路点的电流为disp(Idb););di

38、sp(0s 时各节点的电压标么值为);disp(Vb);disp(0s 时各支路的电流标么值为);disp(Ib)for a=1:NIs=zeros(1,nl);Isy=zeros(1,nl);Idzs=0;s(a)=input( 请输入短路后几秒（若要终止本程序请输入if s(a)=-1disp(本程序已被终止，谢谢使用！);break-1）s=);endfor k=1:nlif B(k,1)=0&B(k,8)=1Sn(k)=B(k,7);cf=Idb./Ib(k);Z1(k)=cf*Z(D(r,1),D(r,1);X(k)=abs(imag(Z1(k);Xjs1(k)=X(k)*Sn(k)

39、./Sb;if s4|Xjs1(k)3.45|B(k,9)=0Is(k)=1./Xjsl(k);if B(k,9)=0Isy(k)=Is(k)*Sb./(sq3*vb(D(r,1); tsl=(无穷大电源 ); ft=num2str(k);ts2=( 对短路点的计算电抗Xjsl，对应的电流标么值，有名值分别为);dn=strcat(tsl.ft,ts2);disp(dn);disp(Xjsl(k);disp(Is(k);disp(Isy(k);elseIsy(k)=Is(k)*Sn(k)./(sq3*vb(D(r,1);tsl=(等值电源 );ft=num2str(k);ts2=( 对短路点的

40、计算电抗Xjsl，对应的电流标么值，有名值分别为);dn=strcat(tsl,ft,ts2);disp(dn);disp(Xjsl(k);disp(Is(k);disp(Isy(k);endelseXjs=Xjs1(k);tsl=(等值电源 );ft=num2str(k);ts2=(对短路点的计算电抗Xjsl 为);dn=strcat(tsl,ft,ts2);disp(dn); disp(Xjs)Is(k)=input( 查与计算电抗Xjs 对应的电流标么值I=);Isy(k)=Is(k)*Sn(k)./(sq3*vb(D(r,1);disp(对应的电流标么值、有名值（单位：KA）分别为 )

41、;disp(Is(k);disp(Isy(k);endIdzs=Idzs+Isy(k);endendId(a)=Idzs;sj=num2str(s(a);ts=(秒对应的短路点短路总电流为（单位： kA）);jg=strcat(sj,ts); disp(jg)disp(Id(a);endplot(s,Id);xlabel( 时间（秒） );ylabel( 短路电流大小的标么值 );title( 短路电流曲线 ); endotherwisebreakend思考题1.计算短路电流的目的是什么？答：选择电气设备,计算各种保护的定值,制定运行方案等2. 在计算短路电流时，是否用到了节点阻抗阵中的所有元

42、素？据此，可以如何简化短路计算程序？答：不需要，可以利用高斯消元法消去部分中间节点和联络节点，化简成等效网络再进行计算。实验四潮流计算1.用牛顿拉夫逊法求解如下非线性方程组，注意保留迭代过程中的中间值。clcclear allsyms x1 x2f1=x1.2+x2.2-5*x1+1;f2=x1.2+x2.2-3*x2-3;J(1,1)=diff(f1,x1);J(1,2)=diff(f1,x2);J(2,1)=diff(f2,x1);J(2,2)=diff(f2,x2);% Xn=input(X0=)Xn=3;3;a=0.01;b=1;t=0;%迭代次数while (b)fn1=subs(f

43、1,x1 x2,Xn(1,1) Xn(2,1); fn2=subs(f2,x1 x2,Xn(1,1) Xn(2,1);fn(1,1)=fn1;fn(2,1)=fn2;Js(1,1)=subs(J(1,1),x1 x2,Xn(1,1) Xn(2,1);Js(1,2)=subs(J(1,2),x1 x2,Xn(1,1) Xn(2,1); Js(2,1)=subs(J(2,1),x1 x2,Xn(1,1) Xn(2,1);Js(2,2)=subs(J(2,2),x1 x2,Xn(1,1) Xn(2,1); Ja=-1*Js;Xns=double(Jafn)Xn=Xns+XnXnsa(1)=abs(X

44、ns(1)Xnsa(2)=abs(Xns(2)Xnsm=max(Xnsa)t=t+1;if Xnsmab=0;else b=1;enddisp(Xn);enddisp(迭代次数 =)disp(t);disp(最终结果为 )disp(Xn);2.上机验证课本例19.1 和例 19.2 的计算结果 ；3.将两种算法的程序计算程序进行合并clcclear;n=5;nl=5;isb=1;pr=0.00001;B1=1 2 0.03i 0 1.05 0;2 3 0.08+0.3i 0.5i 1 0;2 4 0.1+0.35i 0 1 0;3 4 0.04+0.25i 0.5i 1 0;3 5 0.015i 0 1.05 1; B2=0 0 1.05 1.05 0 1;0 3.7+1.3i 1 0 0 2;0 2+1i 1 0 0 2;0 1.6+0.8i 1 0 0 2;5 0 1.05 1.05 0 3;X=1 0;2 0;3 0;4 0;5 0;Y=zeros(n);e=zeros(1,n);f=zeros(1,n);V=zeros(1,n);O=zeros(1,