哈工大电力系统分析短路及潮流计算实验上机程序汇总.docx
《哈工大电力系统分析短路及潮流计算实验上机程序汇总.docx》由会员分享，可在线阅读，更多相关《哈工大电力系统分析短路及潮流计算实验上机程序汇总.docx（17页珍藏版）》请在冰豆网上搜索。
哈工大电力系统分析短路及潮流计算实验上机程序汇总
上机实验
实验一节点导纳阵生成和短路电流计算实验二简单系统的牛顿法潮流计算
姓名：
班级：
141班
学号：
实验说明：
本程序的电路结构来自翁增银、何仰赞主编的《电力系统分析》的例题
实验一节点导纳阵生成和短路电流计算
一、实验目的
根据所给的电力系统，编制短路电流计算程序，通过计算机进行调试，最后完成一个切实可行的电力系统计算应用程序。
通过自己设计电力系统计算程序加深对电力系统分析的理解，同时加强计算机实际应用能力。
二、实验内容
1、编写数据输入、形成导纳阵程序
2、电力系统短路计算实用公式的计算程序及编制和调试
三、实验过程
1、程序代码及说明
%清屏
clear
%读数据
fid=fopen（'node5.txt','r'）;
A=fscanf（fid,'%d',2）;
B=fscanf（fid,'%f',[5,A（2,1）]）;
fclose（fid）;
B=B';
%求节点导纳矩阵
Y=zeros（A（1,1））;%形成5×5的0阵
%%%%%====================================================%导纳阵元素计算
fora=1:
1:
A（2,1）
m=B（a,1）;
n=B（a,2）;
ifB（a,5）>0
Y（m,m）=Y（m,m）+1/（B（a,3）+j*B（a,4））;
Y（n,n）=Y（n,n）+1/（B（a,3）+j*B（a,4））/（B（a,5）^2）;
Y（m,n）=-1/（B（a,3）+j*B（a,4））/B（a,5）;
Y（n,m）=Y（m,n）;
elseifB（a,5）<0
Y（m,m）=Y（m,m）+1/（B（a,3）+j*B（a,4））-j*B（a,5）;
Y（n,n）=Y（n,n）+1/（B（a,3）+j*B（a,4））-j*B（a,5）;
Y（m,n）=-1/（B（a,3）+j*B（a,4））;
Y（n,m）=Y（m,n）;
end
end
%%%%%=======================================================
Y（2,2）=Y（2,2）+j*0.01;
Z=inv（Y）;
If=1/Z（3,3）;%接金属短路，Zf=0
fori=1:
1:
5
V（1,i）=1-Z（i,3）/Z（3,3）;
end
%%%========================================================%计算各支路电流
fora=1:
1:
A（2,1）
m=B（a,1）;
n=B（a,2）;
ifB（a,5）>0
I（m,n）=（B（a,5）*V（1,m）-V（1,n））/（B（a,3）+j*B（a,4））;
elseifB（a,5）<0
I（m,n）=（V（1,m）-V（1,n））/（B（a,3）+j*B（a,4））;
end
end
%%%=======================================================
fid=fopen（'ans.txt','w'）;
fprintf（fid,'Y_matrix\n'）;
%%%=======================================================%输出导纳阵
fori=1:
1:
5
forj=1:
1:
5
k=Y（i,j）;
re=real（k）;
fprintf（fid,'%8.4f',re）;
im=imag（k）;
ifim<0
fprintf（fid,'%8.4fj\t',im）;
elseifim>=0
fprintf（fid,'+%8.4fj\t',im）;
end
end
fprintf（fid,'\n'）;
end
%%%==========================================================
%%%==========================================================%输出阻抗阵，导纳阵求逆
fprintf（fid,'Z_matrix\n'）;
fori=1:
1:
5
forj=1:
1:
5
k=Z（i,j）;
re=real（k）;
fprintf（fid,'%8.4f',re）;
im=imag（k）;
ifim<0
fprintf（fid,'%8.4fj\t',im）;
elseifim>=0
fprintf（fid,'+%8.4fj\t',im）;
end
end
fprintf（fid,'\n'）;
end
%%%=========================================================
%%%=========================================================%打印输出短路电流
fprintf（fid,'If='）;
re=real（If）;
fprintf（fid,'%8.4f',re）;
im=imag（If）;
ifim<0
fprintf（fid,'%8.4fj\t',im）;
elseifim>=0
fprintf（fid,'+%8.4fj\t',im）;
end
fprintf（fid,'\nV\n'）;
%%%========================================================
%%%========================================================%打印输出节点电压
fori=1:
1:
5
fprintf（fid,'Note%d:
V%d=',i,i）;
k=V（1,i）;
re=real（k）;
fprintf（fid,'%8.4f',re）;
im=imag（k）;
ifim<0
fprintf（fid,'%8.4fj\t',im）;
elseifim>=0
fprintf（fid,'+%8.4fj\t',im）;
end
fprintf（fid,'\n'）;
end
%%%========================================================
%%%========================================================%输出打印支路电流
fprintf（fid,'Ibr\n'）;
fori=1:
1:
5
forj=1:
1:
4
k=I（i,j）;
re=real（k）;
im=imag（k）;
if（re~=0||im~=0）
fprintf（fid,'Branch%d--%d:
I%d%d=',i,j,i,j）;
fprintf（fid,'%8.4f',re）;
ifim<0
fprintf（fid,'%8.4fj\n',im）;
elseifim>=0
fprintf（fid,'+%8.4fj\n',im）;
end
end
end
end
%%%==========================================================
fclose（fid）;%关闭文件
附：
'node5.txt'
2、程序输出结果
节点导纳阵：
节点阻抗阵：
短路电流：
If=0.0001+0.1082j
节点电压：
各支路电流：
四、实验总结
这是我的第一次上机实验，感觉稍微有点难，主要还是在工具软件—C语言或者MATLAB的运用上，但是我相信，以后学习中，我会努力掌握的，这是我把理论应用于实际中的必要桥梁！
实验二简单系统的牛顿法潮流计算
一、实验目的
根据所给的电力系统，编制牛顿法潮流计算程序，通过计算机进行调试，最后完成一个切实可行的电力系统计算应用程序。
通过自己设计电力系统计算程序加深对电力系统分析的理解，同时加强计算机实际应用能力。
二、实验内容
电力系统潮流计算的计算程序设计及编制和调试。
三、程序框图
四、实验过程
1、实验程序及说明
clear
fid=fopen（'node4.txt','r'）;%打开输入数据
A=fscanf（fid,'%f',8）;%读8个数
B=fscanf（fid,'%f',[5,A（2,1）]）;
C=fscanf（fid,'%f',[3,（A（1,1）-1）]）;
fclose（fid）;
B=B';
C=C';
B（2,5）=1/B（2,5）;
Y=zeros（A（1,1））;%得4×4的0阵
%%%===============================================================%与前一实验同法求导纳阵
fora=1:
1:
A（2,1）
m=B（a,1）;
n=B（a,2）;
ifB（a,5）>0
Y（m,m）=Y（m,m）+1/（B（a,3）+j*B（a,4））;
Y（n,n）=Y（n,n）+1/（B（a,3）+j*B（a,4））/（B（a,5）^2）;
Y（m,n）=-1/（B（a,3）+j*B（a,4））/B（a,5）;
Y（n,m）=Y（m,n）;
elseifB（a,5）<0
Y（m,m）=Y（m,m）+1/（B（a,3）+j*B（a,4））-j*B（a,5）;
Y（n,n）=Y（n,n）+1/（B（a,3）+j*B（a,4））-j*B（a,5）;
Y（m,n）=-1/（B（a,3）+j*B（a,4））;
Y（n,m）=Y（m,n）;
end
end
end
%%%=================================================================
%%%=================================================================%节点电压赋初值
Q=zeros（2*A（1,1）,1）;
fori=1:
1:
（A（1,1）-1）
ifC（i,3）>0
Q（2*i-1,1）=C（i,3）;
Q（2*i,1）=0;
else
Q（2*i-1,1）=1;
Q（2*i,1）=0;
end
Q（2*A（3,1）-1,1）=A（4,1）;
Q（2*A（3,1）,1）=0;
end
fid=fopen（'answer.txt','w'）;
fprintf（fid,'====================节点电压V=============================\n'）;fprintf（fid,'迭代计数\t\tV1=e1+jf1\t\t\t\tV2=e2+jf2\t\t\t\tV3=e3+jf3\t\t\n'）;
%%%=======================求W阵===================================W=zeros（2*（A（1,1）-1）,50）;
forx=1:
1:
50%设置迭代次数为50次
fori=1:
1:
A（7,1）
k=C（i,1）;
p=0;
q=0;
m=0;
n=0;
forj=1:
1:
A（1,1）
g=real（Y（k,j））;
b=imag（Y（k,j））;
e=Q（2*j-1,1）;
f=Q（2*j,1）;
p=p+g*e-b*f;
q=q+g*f+b*e;
m=m+g*e-b*f;
n=n+g*f+b*e;
end
W（2*k,x）=C（k,2）-Q（2*k-1,1）*p-Q（2*k）*q;
W（2*k-1,x）=C（k,3）-Q（2*k）*m+Q（2*k-1,1）*n;
End
%%%=========================PQ节点=================================
%%%=========================PV节点=================================forl=1:
1:
A（8,1）
k=C（l+A（7,1）,1）;
p=0;
q=0;
m=0;
n=0;
forj=1:
1:
A（1,1）
g=real（Y（k,j））;
b=imag（Y（k,j））;
e=Q（2*j-1,1）;
f=Q（2*j,1）;
p=p+g*e-b*f;
q=q+g*f+b*e;
end
W（2*k,x）=C（k,2）-Q（2*k-1,1）*p-Q（2*k）*q;
W（2*k-1,x）=（C（k,3））^2-（Q（2*k-1,1）^2+Q（2*k）^2）;
End
%%%================================================================
%%%===============================================================%比较是否符合条件
Max=0;
fori=1:
1:
2*（A（1,1）-1）
Max=max（abs（W（i,x））,Max）;
end
ifMaxbreak;
end
%%%=================================================================
%%%=================================================================%求雅克比行列式
fori=1:
1:
A（1,1）-1
%===========================================
forj=1:
1:
A（1,1）-1
%===========================i~=j
ifi~=j
k=C（i,1）;
g=real（Y（k,j））;
b=imag（Y（k,j））;
e=Q（2*k-1,1）;
f=Q（2*k,1）;
J（2*k-1,2*j-1）=-（g*e+b*f）;
J（2*k-1,2*j）=b*e-g*f;
ifi~=3
J（2*k,2*j）=-J（2*k-1,2*j-1）;
J（2*k,2*j-1）=J（2*k-1,2*j）;
else
J（2*k,2*j-1）=0;
J（2*k,2*j）=0;
end
end
%===============================i=j==================ifi==j
k=C（i,1）;
p=0;
q=0;
m=0;
n=0;
g=real（Y（k,k））;
b=imag（Y（k,k））;
J（2*k-1,2*j-1）=-g*Q（2*i-1,1）-b*Q（2*i,1）;
J（2*k-1,2*j）=b*Q（2*i-1,1）-g*Q（2*i,1）;
ifk~=3
J（2*k,2*j-1）=b*Q（2*i-1,1）-g*Q（2*i,1）;
J（2*k,2*j）=g*Q（2*i-1,1）+b*Q（2*i,1）;
end
forr=1:
1:
A（1,1）
g=real（Y（k,r））;
b=imag（Y（k,r））;
e=Q（2*r-1,1）;
f=Q（2*r,1）;
p=p+g*e-b*f;
q=q+g*f+b*e;
end
J（2*k-1,2*j-1）=-p+J（2*k-1,2*j-1）;
J（2*k-1,2*j）=-q+J（2*k-1,2*j）;
ifi~=3
J（2*k,2*j-1）=q+J（2*k,2*j-1）;
J（2*k,2*j）=-p+J（2*k,2*j）;
else
J（2*k,2*j-1）=-2*Q（2*i-1,1）;
J（2*k,2*j）=-2*Q（2*i,1）;
end
end
end
end
%%%==============================================================
%%%==============================================================%解修正方程得修正量Q0
fori=1:
2:
5
J0（i+1,:
）=J（i,:
）;
J0（i,:
）=J（i+1,:
）;
end
J=J0;
Q0=-inv（J）*W（:
x）;
fori=1:
1:
6
Q（i,1）=Q（i,1）+Q0（i,1）;
End
fprintf（fid,'\t%d\t\t',x）;
fori=1:
A（1,1）-1
k=C（i,1）;
ifQ（2*k,1）>=0
fprintf（fid,'%8.4f+%8.4fj\t\t',Q（2*k-1,1）,Q（2*k,1））;
else
fprintf（fid,'%8.4f%8.4fj\t\t',Q（2*k-1,1）,Q（2*k,1））;
end
end
fprintf（fid,'\n'）;
end
%%%=================================================================
%%%==================================平衡点功率======================
k=A（3,1）;
v=0;
j=sqrt（-1）;
forb=1:
A（1,1）
m=conj（Y（k,b））;
p=Q（2*b-1,1）+Q（2*b,1）*j;
n=conj（p）;
v=v+m*n;
end
Wp=（Q（2*k-1,1）+j*Q（2*k,1））*v;
end
%%%===========================节点电压V===========================%fid=fopen（'answer.txt','w'）;
fprintf（fid,'节点电压V\n'）;
fori=1:
A（1,1）-1
k=C（i,1）;
ifQ（2*k,1）>=0
fprintf（fid,'node%d:
%8.4f+%8.4fj\n',k,Q（2*k-1,1）,Q（2*k,1））;
else
fprintf（fid,'node%d:
%8.4f%8.4fj\n',k,Q（2*k-1,1）,Q（2*k,1））;
end
end
fprintf（fid,'平衡点功率P+jQ=%8.4f+%8.4fj',real（Wp）,imag（Wp））;
fclose（fid）;
附注：
node4.txt%输入数据
2、程序输出结果
五、实验总结
本次实验是我把应用理论知识的重要实践，经过实验，我有两点感想，首先，作为一名工科学生，应该能熟练运用C语言和MATLAB等工具，其次，理论如果不用于实践，就永远不知道理论是用来干什么的，学到头一直是满脑子的浆糊，所以，以我的切身经历建议，把这门实验放在跟课程平行的时间上进行，这样不仅有利于实验开展，也有利于学生更加深刻地学习！
/p>
（
关闭
）
2.
在图标菜单点击
全选
和
自动调节缩放3.在功能目录表（图9的X）选择镜像单元4.在形式选择栏选择‘复制’5.在镜像平面选择栏选择‘z-x平面’6.确认y:
‘0’后点击键7.在IconMenu点击全选8.在镜像平面选择栏选择‘y-z平面’9.确认x:
‘0’点击10.点击键ToggleonX图3.9完成模型的腹板部分
开口部详细分析
把临时输入的梁单元扩展为板单元来完成如图11的圆形开口部的翼缘、垂直和水平加劲板以及梁的翼缘。
1.在图标菜单点击
2.图标菜单点击标准视图设定全局坐标系
3.选择工作选项（图10的①）
4.在特性值>截面双击截面号‘999’（曲管形补强加劲板）
5.在单元工具条点击将单元扩展
6.在扩展类型选择栏选择‘线单元→平面单元’
7.在厚度选择栏选择‘10:
0.010000’
8.在生成形式选择栏选择‘复制和移动’
9.在等间距的dx,dy,dz输入栏输入‘0,-0.1,0’
10.在复制次数输入栏输入‘3’输入
11.点击键
选择属性12.在图标菜单点击
13.在选择属性选择‘截面’
14.选择截面号‘998’（垂直、水平加劲板）
15.点击
16.点击键键
17.在厚度选择栏选择‘15:
0.015000’
18.点击键
19.与步骤12∼16相同的方法选择截面号‘997’（梁的翼缘）
20.在厚度选择栏选择‘40:
0.040000’
21.在复制数量输入栏输入‘4’
22.点击
23.点击键键
15
例题
16
升级会员 