已知电力系统中节点1为平衡节点，节点4为PV节点，节点2、3为PQ节点，变压器两侧的电压等级分别为 10kV和 110kV。已知:Ù=1.0520°，P=0.5，U=1.1， Zz=0.08+j0.4，Z3=0.12+j0.5，Zz3=0.1+j0.4，Z=j0.3。S;=100MVA,V=V，负荷及支路阻抗S2=0.6+j0.25,S3=0.25+j0.1,编写详细的matlab潮流计算程序。

以下是一个简单的matlab潮流计算程序，根据输入的参数计算每个节点的电压和功率流量。其中，变量Ybus是节点导纳矩阵，PQ节点和PV节点的注入功率和电压是已知的，通过迭代求解未知节点的电压和功率流量。

% 输入电力系统参数
Zz = 0.08 + j*0.4;
Z3 = 0.12 + j*0.5;
Zz3 = 0.1 + j*0.4;
Z = j*0.3;
Sbase = 100e6;
Vbase = 110e3;
V1 = Vbase;
V4 = 10e3;
theta1 = 0;
theta4 = 0;
P4 = 0.5*Sbase;
U2 = 1.1;
U3 = 1.1;
S2 = 0.6 + j*0.25;
S3 = 0.25 + j*0.1;

% 计算节点导纳矩阵
Y11 = 1/(Zz + Z/2) + 1/(Z3 + Zz3/2);
Y12 = -1/(Zz + Z/2);
Y13 = -1/(Z3 + Zz3/2);
Y14 = 0;
Y21 = -1/(Zz + Z/2);
Y22 = 1/(Zz + Z/2) + 1/Z;
Y23 = 0;
Y24 = -1/Z;
Y31 = -1/(Z3 + Zz3/2);
Y32 = 0;
Y33 = 1/(Z3 + Zz3/2) + 1/Z;
Y34 = -1/Z;
Y41 = 0;
Y42 = -1/Z;
Y43 = -1/Z;
Y44 = 1/Z;

Ybus = [Y11 Y12 Y13 Y14;
        Y21 Y22 Y23 Y24;
        Y31 Y32 Y33 Y34;
        Y41 Y42 Y43 Y44];

% 初始化迭代变量
V = [V1; U2; U3; V4];
P = [0; -real(S2)/Sbase; -real(S3)/Sbase; P4/Sbase];
Q = [0; -imag(S2)/Sbase; -imag(S3)/Sbase; 0];
iter = 0;
max_iter = 100;
tolerance = 1e-6;

% 迭代求解节点电压和功率流量
while iter < max_iter
    iter = iter + 1;

    % 计算注入电流
    I = conj(Ybus*V);

    % 更新节点电压
    for n = 2:3
        V(n) = U(n)*exp(1j*theta(n));
    end

    % 更新未知节点的注入功率和电压
    P(1) = real(V(1)*conj(I(1)));
    Q(1) = imag(V(1)*conj(I(1)));
    U(4) = abs(V(4));
    theta(4) = angle(V(4));
    P(4) = real(V(4)*conj(I(4)));
    Q(4) = imag(V(4)*conj(I(4)));

    % 判断收敛条件
    max_error = max(abs([P; Q]));
    if max_error < tolerance
        break
    end
end

% 输出结果
fprintf('节点电压和功率流量:\n');
fprintf('V1 = %.4f < %.4f\n', abs(V(1)), rad2deg(angle(V(1))));
fprintf('V2 = %.4f < %.4f\n', abs(V(2)), rad2deg(angle(V(2))));
fprintf('V3 = %.4f < %.4f\n', abs(V(3)), rad2deg(angle(V(3))));
fprintf('V4 = %.4f < %.4f\n', abs(V(4)), rad2deg(angle(V(4))));
fprintf('P1 = %.4f MW\n', P(1)*Sbase/1e6);
fprintf('P4 = %.4f MW\n', P(4)*Sbase/1e6);
fprintf('Q1 = %.4f MVar\n', Q(1)*Sbase/1e6);
fprintf('Q4 = %.4f MVar\n', Q(4)*Sbase/1e6);