https://kb.cloudpss.net/documents/software/emtlab/power-flow/fundamentals/

各个部分的参数：

基准容量 (baseMVA)

mpc.baseMVA = 100，表示系统的基准容量为100MVA

母线数据 (bus data) 每行包含以下参数：

bus_i: 母线编号
type: 节点类型（1=PQ节点，2=PV节点，3=平衡节点/参考节点）
Pd: 有功负荷（MW）
Qd: 无功负荷（MVar）
Gs: 对地电导
Bs: 对地电纳
area: 区域号
Vm: 电压幅值初值（标幺值）
Va: 电压相角初值（度）
baseKV: 基准电压（kV）
zone: 区域
Vmax/Vmin: 电压上下限

发电机数据 (generator data)

bus: 发电机连接的母线号
Pg: 有功出力（MW）
Qg: 无功出力（MVar）
Qmax/Qmin: 无功出力上下限
Vg: 发电机端电压设定值
mBase: 发电机容量基准
status: 发电机运行状态（1=运行，0=停运）
Pmax/Pmin: 有功出力上下限

支路数据 (branch data)

fbus/tbus: 支路起始和终止节点
r: 支路电阻（标幺值）
x: 支路电抗（标幺值）
b: 支路对地电纳（标幺值）
rateA/B/C: 支路容量限值（不同时期）
ratio: 变压器变比
angle: 移相角度
status: 支路运行状态（1=运行，0=停运）
angmin/angmax: 相角差限值

matlab


%% 模块1：读取原始数据（ReadData）
[baseMVA, bus, gen, branch, NUM, PL, U, Theta] = ReadData(fileName);


%% 运行Matpower库潮流计算
[results, Vm, Va, Pij_AC, Qij_AC, Sbus, Pbus, Qbus] = runMatpowerPF(fileName);


%% 我的牛顿拉夫逊法计算潮流
tic;  %开始计时
%% 模块2：形成导纳矩阵（ FormYmatrix ）
[Y,G,B] = FormYmatrix(baseMVA,bus,branch,NUM);% 形成导纳矩阵 

%% 牛顿-拉夫逊迭代法求解潮流
% 计算初始节点注入功率
[P,Q] = PQ_inject(baseMVA,bus,gen,NUM);  % 计算各节点的有功、无功功率注入值

% 设置迭代参数
Iteration_Num = 0;           % 当前迭代次数计数器
Iteration_Num_max = 50;      % 最大允许迭代次数，防止死循环
Accuracy = 1e-5;             % 收敛精度要求，即功率不平衡量的允许误差

% 开始牛顿-拉夫逊迭代过程
while Iteration_Num <= Iteration_Num_max    
    % 计算功率不平衡量
    [cort,dP,dQ,Pi,Qi] = CalculateDeltaPQ(P,Q,U,Theta,Y,PL); 
    % cort: 功率不平衡量向量
    % dP: 有功功率不平衡量
    % dQ: 无功功率不平衡量
    % Pi,Qi: 当前迭代步计算的注入功率
    
    max_mismatch = cort;  % 记录当前迭代步的最大功率不平衡量
    
    % 判断是否满足收敛条件
    if ( max(abs(cort)) < Accuracy )  
        fprintf('计算收敛\n');
        fprintf('迭代次数：%d 次\n',Iteration_Num);
        break;  % 满足收敛条件，退出迭代循环
    else
        % 未收敛，继续迭代计算
        % 形成雅可比矩阵
        J = FormJacobi(U,Theta,Y,PL);             
        
        % 求解修正方程，更新状态变量
        [U,Theta] = SolveCorrectionEquation(NUM,U,Theta,PL,J,cort);   
        % U: 更新后的节点电压幅值
        % Theta: 更新后的节点电压相角
        
        % 更新迭代次数
        Iteration_Num = Iteration_Num + 1;                                  
    end
end

% 判断是否因超过最大迭代次数而退出循环
if Iteration_Num > Iteration_Num_max
    fprintf('\n潮流计算不收敛\n');
end

time = toc; %结束计时并记录
%运行时间
fprintf('运行时间：%.4f秒',time);

我来详细解释这个潮流计算程序的每个步骤：

模块1：读取原始数据

matlab
[baseMVA, bus, gen, branch, NUM, PL, U, Theta] = ReadData(fileName);

读取系统基本数据，包括基准容量、母线、发电机、支路信息
初始化电压幅值(U)和相角(Theta)
NUM为系统节点数，PL为节点类型标志

Matpower库潮流计算

matlab
[results, Vm, Va, Pij_AC, Qij_AC, Sbus, Pbus, Qbus] = runMatpowerPF(fileName);

使用Matpower工具箱进行对比验证

自编牛顿拉夫逊法计算潮流

a) 形成导纳矩阵

matlab
[Y,G,B] = FormYmatrix(baseMVA,bus,branch,NUM);

Y是节点导纳矩阵
G是导纳实部（电导）
B是导纳虚部（电纳）

b) 计算初始节点注入功率

matlab
[P,Q] = PQ_inject(baseMVA,bus,gen,NUM);

计算各节点的有功(P)和无功(Q)功率注入值

c) 迭代计算过程：

计算功率不平衡量

matlab
[cort,dP,dQ,Pi,Qi] = CalculateDeltaPQ(P,Q,U,Theta,Y,PL);

检查收敛条件
- 如果最大功率不平衡量小于设定精度，则收敛
- 否则继续迭代
形成雅可比矩阵

matlab
J = FormJacobi(U,Theta,Y,PL);

求解修正方程，更新状态变量

matlab
[U,Theta] = SolveCorrectionEquation(NUM,U,Theta,PL,J,cort);

迭代过程直到：

达到收敛条件
或超过最大迭代次数（50次）

牛顿-拉夫逊法求解电力潮流的理论过程

潮流计算的求解是基于电路中的节点电压方程：

\dot{\mathbf{I}}_n = \mathbf{Y}_n \dot{\mathbf{U}}_n

其中，节点导纳矩阵 $\mathbf{Y}_n$ 描述了网络中元件特性与网络拓扑的连接关系，而节点电压方程则描述了各个节点之间的电流与电压关系。

潮流计算的边界条件是复功率，因此需要建立功率与电压之间的关系，即功率方程（潮流方程）。可以将节点电流用节点功率和电压表示：

\dot{I}_i = \frac{S_i^*}{U_i^*} = \frac{P_i - jQ_i}{U_i^*}

将电流表达式带入节点电压方程后，可以得到潮流方程：

P_i - jQ_i = U_i^* \sum_{j=1}^n Y_{ij} \dot{U}_j \quad (i = 1, 2, \cdots, n)

其中， $Y_{ij} = G_{ij} + jB_{ij}$ ，并且

\dot{U}_i = U_i \angle \delta_i = U_i e^{j\delta_i}

将其带入潮流方程后，可以将潮流方程转化为极坐标形式：

P_i - jQ_i = U_i e^{-j\delta_i} \sum_{j=1}^n \left( G_{ij} + jB_{ij} \right) U_j e^{j\delta_j}

= U_i \sum_{j=1}^n U_j \left( G_{ij} + jB_{ij} \right) e^{-j\delta_{ij}}

= U_i \sum_{j=1}^n U_j \left( G_{ij} + jB_{ij} \right) (\cos\delta_{ij} + j\sin\delta_{ij})

将实部和虚部分开，得到潮流计算的求解方程：

\begin{cases} P_i = U_i \sum_{j=1}^n U_j \left( G_{ij} \cos\delta_{ij} + B_{ij} \sin\delta_{ij} \right) \\ Q_i = U_i \sum_{j=1}^n U_j \left( G_{ij} \sin\delta_{ij} - B_{ij} \cos\delta_{ij} \right) \end{cases}

从以上分析可以看出，一个节点包含 4 个运行变量（ $P$ , $Q$ , $U$ , $\delta$ ），对应 2 个方程。对于 $n$ 个节点的系统，总共有 4n 个变量和 2n 个方程。如果每个节点给定 2 个变量值，那么可以通过这些方程求解剩余的 2n 个未知量。

电力系统节点根据不同的变量给定方式可以分为PQ节点、PV节点和平衡节点。

PQ节点：给定有功功率 $P_i$ 和无功功率 $Q_i$ ，待求节点电压幅值 $U_i$ 和电压相位 $\delta_i$ 。
PV节点：给定有功功率 $P_i$ 和电压幅值 $U_i$ ，待求无功功率 $Q_i$ 和电压相位 $\delta_i$ 。
平衡节点：全系统功率必须平衡，功率损耗 $P_{\text{loss}}$ 和 $Q_{\text{loss}}$ 是状态变量的函数，事先未知。需要一个节点不可以给定 $P$ 和 $Q$ ，但给定电压幅值 $U$ 和相位 $\delta$ ，用以平衡全系统功率，故称为平衡节点。

由于潮流计算方程是一个非线性方程组，直接求解非常困难，通常采用牛顿-拉夫逊法来求解。

牛顿-拉夫逊法求解潮流基本步骤：

形成节点导纳矩阵；
给各节点电压相角赋初值；
代入修正方程，求常数项向量；
求雅可比矩阵；
解修正方程，求解节点电压相角的修正量；
更新节点电压；
检查收敛性，如果不收敛，则用新的节点电压作为初值进行下一次迭代；否则进入下一步；
计算支路功率分布，包括计算PV节点的无功功率和平衡节点的注入功率。

功率不平衡计算：

潮流计算中的功率不平衡列向量为：

\Delta S = \begin{bmatrix} \Delta P \\ \Delta Q \end{bmatrix} = \begin{bmatrix} \Delta P_1 & \cdots & \Delta P_{n-m-1} & \Delta Q_1 & \cdots & \Delta Q_m \end{bmatrix}^T

其中， $n$ 是节点数， $m$ 是PQ节点的数量， $n-m-1$ 是PV节点的数量，1个平衡节点。

对于PQ节点的功率不平衡量：

\Delta P_i = P_i - |V_i| \sum_{j=1}^n |V_j| \left( G_{ij} \cos{\theta_{ij}} + B_{ij} \sin{\theta_{ij}} \right)

\Delta Q_i = Q_i - |V_i| \sum_{j=1}^n |V_j| \left( G_{ij} \sin{\theta_{ij}} - B_{ij} \cos{\theta_{ij}} \right)

其中， $i = 1, 2, \dots, m$ ， $V_i$ 是电压幅值， $\theta_{ij}$ 是电压相位差。

对于PV节点：

电压幅值已给定，因此不作为变量。在迭代过程中，只考虑PV节点的有功功率偏差方程：

\Delta P_i = P_i - |V_i| \sum_{j=1}^n |V_j| \left( G_{ij} \cos{\theta_{ij}} + B_{ij} \sin{\theta_{ij}} \right)

其中， $i = 1, 2, \dots, n-m-1$ 。

雅可比矩阵的计算：

雅可比矩阵 $J$ 的元素可以分为4个部分，分别是 $H$ 、 $N$ 、 $J$ 和 $L$ ，具体计算如下：

非对角元素：

H_{ij} = -V_i V_j \left( G_{ij} \sin{\theta_{ij}} - B_{ij} \sin{\theta_{ij}} \right)

N_{ij} = -V_i V_j \left( G_{ij} \cos{\theta_{ij}} + B_{ij} \sin{\theta_{ij}} \right)

J_{ij} = V_i V_j \left( G_{ij} \cos{\theta_{ij}} + B_{ij} \sin{\theta_{ij}} \right)

L_{ij} = -V_i V_j \left( G_{ij} \sin{\theta_{ij}} - B_{ij} \cos{\theta_{ij}} \right)

对角元素：

H_{ii} = Q_i + B_{ii} V_i^2

N_{ii} = -P_i - G_{ii} V_i^2

J_{ii} = -P_i + G_{ii} V_i^2

L_{ii} = -Q_i + B_{ii} V_i^2

牛顿-拉夫逊法迭代计算：

在牛顿-拉夫逊法中，通过不断迭代计算，直到所有节点的修正量满足收敛条件。

修正方程为：

\begin{bmatrix} \Delta P \\ \Delta Q \end{bmatrix} = \begin{bmatrix} H & N \\ J & L \end{bmatrix} \begin{bmatrix} \Delta \theta \\ \Delta V^{-1} \Delta V \end{bmatrix}

通过解得修正量 $\Delta \theta$ 和 $\Delta V$ ，更新节点电压：

\begin{bmatrix} \theta^{(k+1)} \\ V^{(k+1)} \end{bmatrix} = \begin{bmatrix} \theta^{(k)} \\ V^{(k)} \end{bmatrix} + \begin{bmatrix} \Delta \theta^{(k+1)} \\ \Delta V^{(k+1)} \end{bmatrix}

反复迭代，直到满足收敛条件： $\Delta V < \varepsilon$ 和 $\Delta \theta < n$ ，其中 $\varepsilon$ 是允许的误差， $n$ 是最大迭代次数。

目录