实验指导书.docx
《实验指导书.docx》由会员分享，可在线阅读，更多相关《实验指导书.docx（20页珍藏版）》请在冰点文库上搜索。
实验指导书
内蒙古工业大学电力学院
《电力系统计算机辅助潮流计算》
实验指导书
编者：
顾秀芳
单位：
电力学院电力工程系
时间：
2006年5月
《电力系统计算机辅助潮流计算》实验指导书
1、实验目的：
了解计算机潮流分析的基本原理、主要步骤；掌握节点导纳矩阵形成和修改的方法，掌握数据处理的基本方法；熟悉Matlab运行环境，了解Matlab基本编程语句和语法；运用潮流分析程序对给定网络的运行方式做潮流分析，并初步分析计算结果；
2、实验内容：
输入网络参数，包括节点号、节点导纳矩阵、节点功率等；输入潮流程序、调试并输出结果，绘制潮流分布图。
3、实验要求：
通过预习，对计算机潮流分析基本理论有深入了解；为程序准备必要的、准确的原始数据；熟悉Matlab运行环境，输入潮流程序，上机独立完成程序的调试，给出潮流分析的结果并按要求绘制潮流分布图。
4、实验步骤：
（1）熟悉原始资料：
根据计算要求，整理数据，包括：
计算网络中线路、变压器的参数、形成节点导纳矩阵；表示各节点的注入功率。
（以上数据均采用有名值计算）
（2）读通潮流程序：
完成程序的解释和说明，必要时附加对应的公式和程序语言的说明；
（3）上机调试：
熟悉Matlab的运行环境，准确输入原始数据、节点编号、节点注入功率等信息；
（4）整理计算结果：
根据计算结果作电网潮流分布图。
5、原始数据：
（1）网络接线示意图；
（2）导线、变压器型号及参数；
（3）节点负荷或发电机功率。
6、实验题目：
1）实验网络接线图及原始数据：
如图所示，3为平衡节点，1、2为P、Q节点，电压等级为110kV，节点处功率已将各线路充电功率考虑在内，3节点电压为115kV，角度为0。
表1-1导线型号参数的选择
线路型号
（
）
（
）
线路长度（
）
长期载流量A
1-2
LGJ-70
0.45
0.441
53.3
297
1-3
LGJ-150
0.21
0.416
47.62
472
2-3
LGJ-150
0.21
0.416
64.29
472
2）要求完成的内容：
（1）预习报告：
按照程序要求，对图中节点重新编号，计算各节点自导纳、互导纳，确定节点注入功率，分析程序中需要修改的地方。
（2）上机实验一：
重点熟悉MATLAB的基本操作，基本语句的含义；
（3）上机实验二：
在MATLAB环境下导入程序，输入原始数据，取允许误差ε=10-4；
（4）上机调试：
完成正常运行的潮流分析
1节点导纳矩阵；
2第一次迭代用的雅可比矩阵，各节点电压的变化量、功率不平衡量（正常潮流）；
3迭代收敛后各节点电压，支路功率，迭代次数，平衡节点功率；
（5）分析事故断线下潮流分布（只列出最后迭代结果），并判断线路的过载情况；
（6）实验报告需要完成：
正常运行及断线故障下的潮流分布图的绘制（潮流的标注请参考图所示）；
潮流分布示意图范例：
并回答以下问题：
①如何检查Y的正确性；
②如何判断潮流结果的正确性与合理性。
7、参考文献：
（1）陈跃主编，《电气工程专业毕业设计指南-电力系统分册》，中国水利水电出版社，2003年；
（2）李敏，代莹，陈金富，段献忠.潮流计算收敛问题研究综述.继电器.2006，11（3）：
74-76
附录一：
潮流计算程序
程序说明：
此潮流程序选自《电气工程专业毕业设计指南-电力系统分册》中潮流计算程序，应用于计算含平衡节点、PQ节点的电力网的潮流计算及分析。
程序要求整理对原始数据进行简单的处理后方可应用。
需要在程序应用前准备：
网络节点导纳矩阵、PQ节点数据等。
原始网络：
Clear%程序数据以上给定网络的数据，实验中需要输入实际的网络参数
G（1,1）=12.900;B（1,1）=-31.400;G（1,2）=-7.300;B（1,2）=22.5;
G（1,3）=-2.400;B（1,3）=4.700;G（1,4）=-3.200;B（1,4）=4.200;
G（1,5）=0;B（1,5）=0;
G（2,1）=-7.300;B（2,1）=22.5;G（2,2）=7.300;B（2,2）=-22.500;G（2,3）=0;B（2,3）=0;
G（2,4）=0;B（2,4）=0;G（2,5）=0;B（2,5）=0;
G（3,1）=-2.400;B（3,1）=4.700;G（3,2）=0;B（3,2）=0;G（3,3）=4.900;B（3,3）=-7.900;
G（3,4）=-2.500;B（3,4）=3.200;G（3,5）=0;B（3,5）=0;
G（4,1）=-3.200;B（4,1）=4.200;G（4,2）=0;B（4,2）=0;G（4,3）=-2.500;B（4,3）=3.200;
G（4,4）=9.4600;B（4,4）=-12.250;G（4,5）=-3.76;B（4,5）=4.8500;
G（5,1）=0;B（5,1）=0;G（5,2）=0;B（5,2）=0;G（5,3）=0;B（5,3）=0;
G（5,4）=-3.760;B（5,4）=4.850;G（5,5）=3.760;B（5,5）=-4.850;%分别求各互导纳、自导纳，若题目给出阻抗，则需要先转换为导纳后计算；输入计算中也可利用导纳的对称性简化输入，如G（2，3）=G（3，2）
Y=G+j*B;%形成节点导纳矩阵；
delt
（1）=0;delt
（2）=0;delt（3）=0;delt（4）=0;
u
（1）=1.0;u
（2）=1.0;u（3）=1.0;u（4）=1.0;%设置P、Q节点的幅值和角度，若采用有名值计算，电压幅值采用额定电压数值，角度不论采用有名值或标幺值通常均设为0；
p
（1）=0.62;q
（1）=0.30;p
（2）=-0.1734;q
（2）=-0.084;
p（3）=-0.2231;q（3）=-0.105;p（4）=-0.15;q（4）=-0.0728;为P、Q节点的注入有功和注入无功赋给定值；
k=0;precision=1;设置迭代次数初始值k，赋精度初始值，此值应大于计算要求的精度，以进入以下循环；
N1=4;%P、Q节点的数量；
whileprecision>0.00001%判断精度是否满足要求，若不满足，则继续迭代
delt（5）=0;u（5）=1.00;%平衡节点编号在程序设为最大，此处将平衡节点的电压幅值和角度赋值；
clear
G（1,1）=10.834;B（1,1）=-32.400;G（1,2）=-1.667;B（1,2）=5;
G（1,3）=-1.667;B（1,3）=5;G（1,4）=-2.5;B（1,4）=7.5;
G（1,5）=-5;B（1,5）=15;
G（2,1）=-1.667;B（2,1）=5;G（2,2）=12.917;B（2,2）=-38.75;G（2,3）=-10;B（2,3）=30;
G（2,4）=0;B（2,4）=0;G（2,5）=-1.25;B（2,5）=3.75;
G（3,1）=-1.667;B（3,1）=5;G（3,2）=-10;B（3,2）=30;G（3,3）=12.917;B（3,3）=-38.75;
G（3,4）=-1.25;B（3,4）=3.75;G（3,5）=0;B（3,5）=0;
G（4,1）=-2.5;B（4,1）=7.5;G（4,2）=0;B（4,2）=0;G（4,3）=-1.25;B（4,3）=3.75;
G（4,4）=3.75;B（4,4）=-11.25;G（4,5）=0;B（4,5）=0;
G（5,1）=-5;B（5,1）=15;G（5,2）=-1.25;B（5,2）=3.75;G（5,3）=0;B（5,3）=0;
G（5,4）=0;B（5,4）=0;G（5,5）=6.25;B（5,5）=-18.75;
Y=G+j*B%形成节点导纳矩阵
delt
（1）=0;delt
（2）=0;delt（3）=0;delt（4）=0;
u
（1）=1.0;u
（2）=1.0;u（3）=1.0;u（4）=1.0;
p
（1）=0.2;q
（1）=0.20;p
（2）=-0.45;q
（2）=-0.15;
p（3）=-0.4;q（3）=-0.05;p（4）=-0.6;q（4）=-0.1;%设迭代初值
k=1;precision=1;%设迭代次数和精度
N1=4;whileprecision>0.0001%判断是否满足精度要求
form=1:
N1
forn=1:
N1+1
pt（n）=u（m）*u（n）*（G（m,n）*cos（delt（m）-delt（n））+B（m,n）*sin（delt（m）-delt（n）））;
qt（n）=u（m）*u（n）*（G（m,n）*sin（delt（m）-delt（n））-B（m,n）*cos（delt（m）-delt（n）））;
end
pi（m）=sum（pt）,qi（m）=sum（qt）;%计算PQ节点的注入功率
dp（m）=ps（m）-pi（m）;
dq（m）=qs（m）-qi（m）;%计算PQ节点的功率不平衡量
end
form=1:
N1
forn=1:
N1+1
pt（n）=u（m）*u（n）*（G（m,n）*cos（delt（m）-delt（n））+B（m,n）*sin（delt（m）-delt（n）））;
%
程序中m
i，n
j，，j=1,2,….N+1,
qt（n）=u（m）*u（n）*（G（m,n）*sin（delt（m）-delt（n））-B（m,n）*cos（delt（m）-delt（n）））;
%
，
end
pi（m）=sum（pt）,qi（m）=sum（qt）;
dp（m）=ps（m）-pi（m）;
dq（m）=qs（m）-qi（m）;
end
form=1:
N1
forn=1:
N1
ifm==n
H（m,m）=-qi（m）-u（m）^2*B（m,m）;N（m,m）=pi（m）+u（m）^2*G（m,m）;
J（m,m）=pi（m）-u（m）^2*G（m,m）;L（m,m）=qi（m）-u（m）^2*B（m,m）;
JJ（2*m-1,2*m-1）=H（m,m）;JJ（2*m-1,2*m）=N（m,m）;
JJ（2*m,2*m-1）=J（m,m）;JJ（2*m,2*m）=L（m,m）;
else
H（m,n）=u（m）*u（n）*（G（m,n）*sin（delt（m）-delt（n））-B（m,n）*cos（delt（m）-delt（n）））;
J（m,n）=-u（m）*u（n）*（G（m,n）*cos（delt（m）-delt（n））+B（m,n）*sin（delt（m）-delt（n）））;
N（m,n）=-J（m,n）;L（m,n）=H（m,n）;
JJ（2*m-1,2*n-1）=H（m,n）;JJ（2*m-1,2*n）=N（m,n）;
JJ（2*m,2*n-1）=J（m,n）;JJ（2*m,2*n）=L（m,n）;
%
end
end
end%计算jocbi各项,并放入统一矩阵JJ中，对JJ下标统一编号
JJ
form=1:
N1
PP（2*m-1）=pp（m）;
PP（2*m）=qq（m）;
end%按统一矩阵形成功率不平衡
uu=inv（JJ）*PP';precision=max（abs（uu））;%判断是否收敛
forn=1:
N1
delt（n）=delt（n）+uu（2*n-1）;
u（n）=u（n）+uu（2*n）*u（n）;%将结果分解为电压幅值和角度
end%求解修正方程，得电压幅值变化量（标幺值）和角度变化量
k=k+1;
end
forn=1:
N1+1
U（n）=u（n）*（cos（delt（n））+j*sin（delt（n）））;
end
form=1:
N1+1
I（m）=Y（5,m）*U（m）;%求平衡节点的注入电流
end
S5=U（5）*sum（conj（I））%求平衡节点的注入功率
form=1:
N1+1
forn=1:
N1+1
S（m,n）=U（m）*（conj（U（m））-conj（U（n）））*conj（-Y（m,n））;
end
end%求节点i,j节点之间的功率，方向为由i指向j,
S%显示支路功率
附录二：
Matlab基本程序语言介绍
一、MATLAB简介
MATLAB是矩阵（Matric）和实验室（Laboratory）两个英文单词的前三个字母的组合。
它的首创者是美国新墨西哥大学计算机系的系主任CleverMoler博士。
现在MATLAB已经被广泛应用于数值计算、图形处理、符号运算、数学建模、小波分析、系统辨识、实时控制和动态仿真等研究领域。
在国外，MATLAB已经成为大学生及研究生所必须掌握的一门计算机语言，并且成为科学研究、工程设计与运算等任务的得力助手。
二、MATLAB语言环境介绍
1、菜单栏区
菜单栏区主要包括文档（File）菜单、编辑（edit）菜单、视图（View）菜单、网络（Window）菜单和帮助（Help）菜单。
2、常用工具栏区
常用工具栏主要包括新建M文档按钮、打开按钮、剪切按钮、复制按钮等。
可以将鼠标箭头放在按钮上识别。
3、窗口区
窗口区主要包括工作区（Workspace）、指令窗口（CommandWindow）和指令历史窗口（CommandHistory）。
4、开始导航区
开始导航区由开始（Start）按钮（在屏幕左下方）进行引导。
利用导航区可以运行所有MATLAB软件的工具以及访问相关文档。
三、MATLAB基本命令及语句
1、MATLAB输入命令的常用方式有两种：
①直接在MATLAB的命令窗口中逐条输入MATLAB命令，适用于命令行很简单时；②m文件工作方式，适用于命令行很多时，易于把握程序的具体走向，易于程序的修改和维护。
每次运行时只需在MATLAB的命令窗口输入m文件的文件名就可以了。
注意应将m文件放在MATLAB的默认搜索路径下（X：
\MATLAB6P1\WORK,X表示MATLAB所在的硬盘）,这样就不用设置文件的路径了,否则应当用路径操作指令path重新设置路径。
注意：
MATLAB中变量、文件名、函数名必须以字母开头，可以由字母、数字和下划线混合组成，名称区分大小写，尽量避免特定意义的默认变量，如：
i,j,pi,ans,inf等。
例1：
在MATLAB桌面执行“File/New/M-file”命令，或单击“New/M-file”快捷键按钮。
例：
A=[123;456;780];
B=[230;587;1036];
C=A+B
invC=inv（C）
输入完毕，单击保存按钮，键入“test”，建立test.m文件（注意路径）
运行方法：
（1）在命令窗口键入“test”；
（2）在输入M-file窗口中，选“Debug-Run”。
观察MATLAB工作空间Workplace变量。
运行结果:
C=
353
91313
17116
invC=
-0.23720.01090.0949
0.6095-0.1204-0.0438
-0.44530.1898-0.0219
数据格式的显示:
A=[1,2;4,5];
B=[3;5];
X=A\B
1）数据格式选择:
short:
显示5位定点格式;long:
显示15位定点格式;shorte:
显示5位浮点格式;long:
显示15位浮点格式
2）数字显示选项:
loose:
语句之间加入间隔行,增加语句的可读性;compact:
压缩多余的行,使得窗口中显示出更多的语句行.
2、矩阵的形成：
两种方法：
①行中各元素无符号：
例：
B=[123；234]
②行中各元素以逗号隔开：
例：
B=[1，2，3；2，3，4]
书写语句时，分号用在每行命令的结尾，则MATLAB执行完该行命令时不会显示本行的计算结果。
（上机试验）当分号用在[]内时，表示矩阵中该行的结束；
若矩阵未指明变量，则在CommandWindow中以ans作为默认变量。
矩阵的使用：
例：
C=[353;91313;17116];
D=C（2，2）%提取第2行第2列的元素
E=C（[2,3],[1:
3]）%提取第2、3行，第1列到第3列的子矩阵；
C（2，：
）=[123]%对C矩阵的第2行赋值，注意赋值的矩阵的维数，与要求赋值的结构相对应。
C（2，：
）=[]%删除第二行的元素
3、常用函数
*复数的形成:
1.B=a+j*b2.B=a+i*b
*复数矩阵的形成:
两个实数矩阵相乘得到:
例:
A=[2,3;4,5],B=[2,4;5,1],C=A+j*B
运行结果:
A=
23
45
B=
24
51
C=
2.0000+2.0000i3.0000+4.0000i
4.0000+5.0000i5.0000+1.0000i
*sum函数：
进行矩阵元素的求和。
格式：
B=sum（A），B=sum（A，dim）
B=sum（A）：
计算每一列元素求和；若是1*N（或N*1）阶矩阵，则将各元素相加。
例1：
C=[353;91313;17116];
D=sum（C）
结果:
C=
353
91313
17116
例2：
A=[123]；B=sum（A）
例3：
A=[1；2；3]；B=sum（A）
*转置：
B=Aˊ
例：
A=[123；456]；B=Aˊ
*实数的绝对值以及复数的模：
B=abs（A）
返回数组或矩阵相应元素的绝对值。
如果是复数矩阵，则返回复数的模值。
若复数为：
A=C+jB，D=abs（A），其中D=（C2+B2）1/2
例1:
实数数组:
A=[1,2,-3];B=abs（A）结果:
B=123
例2:
复数矩阵:
A=[1+j*2,2-j*3;1-j*2,-3-j*4];B=abs（A）
运行结果:
B=
2.23613.6056
2.23615.0000
*求复数的共轭：
B=conj（A）
例1：
A=[1+j2]
B=conj（A）
例2:
A=[1+j*2,1-j*2;1-j*3,1-j*4];B=conj（A）
运行结果:
B=
1.0000-2.0000i1.0000+2.0000i
1.0000+3.0000i1.0000+4.0000i
*矩阵最大元素和最小元素
含义：
在MATLAB软件中，利用max函数求解矩阵的最大元素，用min求解最小元素。
以max为例说明：
格式：
C=max（A），C=max（A，B）
其中：
C=max（A）函数用来求解单一矩阵A的最大元素。
例1：
A=[1，2，3，4]
C=max（A）
例2：
A=[1；2；3；4]
C=max（A）
例1、2结果相同。
C=max（A，B）函数用来求解矩阵A、B中的最大元素。
例：
A=[1，9；8，2]
B=[3，4；6，5]
C=max（A，B）
运行结果:
c=
39
85
*+、—、*：
矩阵的加、减、乘（注意：
+、—运算需要满足矩阵具有相同的维数；A*B乘法运算能够运算的前提条件是A的列数与B的行数相同）
*矩阵的除法：
矩阵右除A/B表示的是方程X*B=A的解A×B-1,矩阵左除A\B表示方程A*X=B的解A-1×B.
例1：
执行以下命令：
A=[1，2；3，4]；
B=[3，5；7，9]；
X=A/B
X是XB=A的解。
例2：
执行以下程序：
A=[1,2;4,5];
B=[3;5];
X=A\B
X是AX=B的解。
运行结果:
X=
-1.6667
2.3333
*求逆：
B=inv（A）,注意求逆的矩阵必须是方阵.
例：
A=[12;45];B=[3;5];X=inv（A）*B
上述求解AX=B的解也可以用实现。
运行结果:
X=
-1.6667
2.3333
与X=A\B相同
*矩阵的形成：
两种方法：
①行中各元素无符号：
例：
B=[123；234]
②行中各元素以逗号隔开：
例：
B=[1，2，3；2，3，4]
基本语句练习题目:
求解方程的解:
x-y+z+2t=1
x+y-2z+t=1
x+y+t=2
x+z-t=1
程序的编写：
a=[1,-1,1,2;1,1,-2,1;1,1,0,1;1,0,1,-1]
b=[1;1;2;1]或：
b=[1,1,2,1]'
c=inv（a）*b
运行结果：
c=
0.8000
0.9000
0.5000
0.3000
*求矩阵的行、列数：
d=size（A）,[m,n]=size（A）
将矩阵A的行数和列数赋给变量d或m，n
*求矩阵的长度：
d=length（A）
计算矩阵的长度（列数）
四、获取帮助命令
1）help命令
·help：
列出最原始的帮助项目，每项对应MATLAB的一个目录名；
·helpfunc：
根据用户指定的项目，func给出相应的帮助信息。
注意：
help必须小写，func查找项目不区分大小写。
2）查询命令
·lookforabc：
系统按照已设置的MATLAB路径在所有文件中查找字符串。
3）who：
列出当前工作空间中的变量；
4）whos：
列出当前工作空间中的变量的详细信息；
5）%：
注释符，在一行的某处键入%，则该行%以后的部分将视维注释部分而不再执行；
6）clear：
清理内存变量；clc:
清除工作窗中的所有显示内容；...：
续行符号。
7）：
：
既可用作矩阵的下标，又可以应用于行循环操作等。
j：
k：
等价于[jj+1j+2…..k]
矩阵的关系运算：
>:
大于;<:
小于;==:
等于;>=大于或等于;<=:
小于或大于;~=:
不等于
3、循环语句：
1）for-end循环
一般形式：
Forvariable=expression
statement
end
2）条件转移结构
if条件式
条件块语句组1
else
条件块语句组2
end
条件语句的另一形式：
if条件式1
条件块语句组1
elseif条件式2
条件块语句组2
…….
elseif
条件块语句组n-1
else
条件块语句组n
end
·while-end循环
其表达式为：
while表达式
循环体语句
end
四、程序的调试
调试程序的步骤：
调试程序的步骤：
1、在程序可疑处设置断点；
2、执行程序；
3、检查程序运行至断点处的各个变量的值；
4、进行处理；
5、恢复执行；
6、结束调试状态；
7、取消断点。
进入建立命令文件的窗口，有Debug菜单和Breakpoint菜单，这两个菜单是专为调试程序提供的。
Breakpoint菜单中选项：
Set/ClearB
升级会员 