end
%----------------- 以上为求雅可比矩阵的各个元素及扩展列的功率差或电压差 ----------------------------
for k=3:N0 % N0=2*n （从第三行开始，第一、二行是平衡节点）
k1=k+1;N1=N; % N=N0+1 即 N=2*n+1扩展列△P、△Q 或 △U
for k2=k1:N1 % 从k+1列的Jacobi元素到扩展列的△P、△Q 或 △U
J(k,k2)=J(k,k2)./J(k,k);% 用K行K列对角元素去除K行K列后的非对角元素进行规格化 end
J(k,k)=1; % 对角元规格化K行K列对角元素赋1 %---------------------- 回代运算 ----------------------------------------
if k~=3 % 不是第三行 k > 3 k4=k-1;
for k3=3:k4 % 用k3行从第三行开始到当前行的前一行k4行消去
for k2=k1:N1 % k3行后各行上三角元素
J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算（当前行k列元素消为0）
end %用当前行K2列元素减去当前行k列元素乘以第k行K2列元素
J(k3,k)=0; %当前行第k列元素已消为0 end
if k==N0 %若已到最后一行 break; end
%---------------------- 前代运算 ---------------------------
for k3=k1:N0 % 从k+1行到2*n最后一行
for k2=k1:N1 % 从k+1列到扩展列消去k+1行后各行下三角元素
J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算
end %用当前行K2列元素减去当前行k列元素乘以第k行K2列元素
J(k3,k)=0; %当前行第k列元素已消为0 end
else %是第三行k=3
%---------------------- 第三行k=3的前代运算 -----------------------------------------
for k3=k1:N0 %从第四行到2n行（最后一行）
for k2=k1:N1 %从第四列到2n+1列（即扩展列）
J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算（当前行3列元素消为0）
end %用当前行K2列元素减去当前行3列元素乘以第三行K2列元素
J(k3,k)=0; %当前行第3列元素已消为0 end end end
%--------上面是用线性变换方式高斯消去法将Jacobi矩阵化成单位矩阵------- for k=3:2:N0-1 L=(k+1)./2;
e(L)=e(L)-J(k,N); %修改节点电压实部 k1=k+1;
f(L)=f(L)-J(k1,N); %修改节点电压虚部 end
%------修改节点电压----------- for k=3:N0
DET=abs(J(k,N));
if DET>=pr %电压偏差量是否满足要求 IT2=IT2+1; %不满足要求的节点数加1 end end
ICT2(a)=IT2; %不满足要求的节点数 ICT1=ICT1+1; %迭代次数 end
%用高斯消去法解\disp('迭代次数：'); disp(ICT1);
disp('没有达到精度要求的个数：'); disp(ICT2); for k=1:n
V(k)=sqrt(e(k)^2+f(k)^2); %计算各节点电压的模值
sida(k)=atan(f(k)./e(k))*180./pi; %计算各节点电压的角度 E(k)=e(k)+f(k)*j; %将各节点电压用复数表示 end
%------------------------ 计算各输出量 ------------------------------- disp('各节点的实际电压标幺值E为(节点号从小到大排列)：'); disp(E); %显示各节点的实际电压标幺值E用复数表示 disp('-----------------------------------------------------');
disp('各节点的电压大小V为(节点号从小到大排列)：'); disp(V); %显示各节点的电压大小V的模值 disp('-----------------------------------------------------');
disp('各节点的电压相角sida为(节点号从小到大排列)：');
disp(sida); %显示各节点的电压相角 for p=1:n C(p)=0; for q=1:n
C(p)=C(p)+conj(Y(p,q))*conj(E(q)); %计算各节点的注入电流的共轭值 end
S(p)=E(p)*C(p); %计算各节点的功率 S = 电压 X 注入电流的共轭值 end
disp('各节点的功率S为(节点号从小到大排列)：'); disp(S); %显示各节点的注入功率
五、实验结果
在上图所示的简单电力系统中，系统中节点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.500172 2. 0.958690 -0.108387 0.964798 -6.450306 3. 1.092415 0.128955 1.100000 6.732347 4. 1.050000 0.000000 1.050000 0.000000
各节点功率：
节点 P Q 1 -0.300000 -0.180000 2 –0.550000 -0.130000 3 0.500000 -0.551305 4 0.367883 0.264698
导纳矩阵 Y=
1.0421 - 8.2429i -0.5882 + 2.3529i 0.0000 + 3.6667i -0.4539 + 1.8911i -0.5882 + 2.3529i 1.0690 - 4.7274i 0.0000 + 0.0000i -0.4808 + 2.4038i 0.0000 + 3.6667i 0.0000 + 0.0000i 0.0000 - 3.3333i 0.0000 + 0.0000i -0.4539 + 1.8911i -0.4808 + 2.4038i 0.0000 + 0.0000i 0.9346 - 4.2616i
迭代次数： 4
各节点的实际电压标幺值E为(节点号从小到大排列)：
0.9854 - 0.0485i 0.9691 - 0.0759i 1.0968 + 0.0844i 1.0500 + 0.0000i 各节点的电压大小e为(节点号从小到大排列)： 0.9854 0.9691 1.0968 1.0500 各节点的电压大小f为(节点号从小到大排列)： -0.0485 -0.0759 0.0844 0 各节点的电压大小V为(节点号从小到大排列)： 0.9866 0.9721 1.1000 1.0500 各节点的电压相角sida为(节点号从小到大排列)： -2.8167 -4.4784 4.4018 0
各节点的功率S为(节点号从小到大排列)：
-0.5500 - 0.1300i -0.3000 - 0.1800i 0.5000 + 0.0857i 0.3595 + 0.2343i
实验六 水火电厂间有功功率负荷的经济分配
一、实验目的
本实验通过对水火电厂间有功功率负荷的经济分配的计算机程序编制与调试，获得能计算有功负荷经济分配的计算机程序，即根据已知的水火电厂的耗量特性、发电用水量及各时间段的有功负荷功率，由计算机程序自动完成各时段各电厂的有功出力的计算。通过实验教学加深学生对电力系统经济运行概念的理解，学会运用数学模型进行水火电厂间的有功功率负荷的经济分配计算，熟悉各种常用应用软件，熟悉硬件设备的使用方法，加强编制调试计算机程序的能力，提高工程计算的能力，学习如何将理论知识和实际工程问题结合起来。
二、实验器材
计算机、软件（已安装，包括各类编程软件C语言、C++、VB、VC等、应用软件MATLAB等）、移动存储设备（学生自备，软盘、U盘等）
三、实验内容
编制调试水、火电厂间有功功率负荷经济分配的计算机程序。程序要求根据已知参数，完成水火电厂间有功功率负荷的经济分配计算，要求计算出各电厂满足经济运行要求所分配的有功负荷等参数。
1.熟悉电力系统经济运行的计算方法，按照计算方法编制程序。 2.将事先编制好的计算程序原代码由自备移动存储设备导入计算机。 3.在相应的编程环境下对程序进行组织调试。 4.应用计算例题验证程序的计算效果。 5.对调试正确的计算程序进行存储、打印。 6.完成本次实验的实验报告。
四、实验数据
一个火电厂和一个水电厂并联运行。火电厂的燃料消耗特性为：
水电厂的耗水量特性为：
水电厂的给定日用水量为
。系统的日负荷变化如下：
时，负荷为700MW，
时 负荷为
时 负荷为350MW，
500MW。火电厂容量为600MW，水电厂容量为450MW。试确定水、火电厂间的功率经济分配。