谁能帮我指教一下改怎么编潮流计算

zhouyou529

clear;
, A) G6 z* B# v4 L; B/ L0 g7 K) B7 w8 [Yb=[6.25-18.75i  -5+15i        -1.25+3.75i    0             0
# o2 q0 c8 S# S) }- k! R    -5+15i       10.834-32.5i  -1.667+5i     -1.667+5i    -2.5+7.5i
1 p) H0 u6 z, D! P% v1 g/ D" q    -1.25+3.75i  -1.667+5i     12.917-38.75i -10+30i      0
% K  r- |# _) j* p0 K. Q    0            -1.667+5i     -10+30i       12.917-38.75i  -1.25+3.75i
    0            -2.5+7.5i     0             -1.25+3.75i    3.75-11.25i];
2 y( F4 V: t, P6 A) KPknow=[0;0.2;-0.45;-0.4;-0.6];Qknow=[0;0.2;-0.15;-0.05;-0.1];e1=1.06;f1=0;%已知量
* R# x  M, k& ?0 T# d! ce=[1.06;1;1;1;1];f=[0;0;0;0;0];
k=0;%跌代次数
PQX=[1;1;1;1;1;1;1;1];%修正量   
! q# ?8 Z. Y) C8 @- w- MP=zeros(4,1);Q=zeros(4,1);I=zeros(4,1);YAKEBI=zeros(8,8);
2 j' d  Z4 k8 O% t1 Qwhile (max(PQX)>=0.1|min(PQX)<=-0.1)
% m2 x4 J% a9 G. y+ |) j* c    k=k+1;
4 l+ U) R4 Y2 c6 z3 q    for n=1:4
( J9 j. }1 G0 A1 l        for m=1:5
/ Q% u: _! o* Q$ s- e( b% s5 [            P(n)=P(n)+e(n+1)*(real(Yb(n+1,m))*e(m)-imag(Yb(n+1,m))*f(m))+f(n+1)*(real(Yb(n+1,m))*f(m)+imag(Yb(n+1,m)*e(m)));
+ J- B) i9 A  Y- K" G" `9 Q            Q(n)=Q(n)+f(n+1)*(real(Yb(n+1,m))*e(m)-imag(Yb(n+1,m))*f(m))-e(n+1)*(real(Yb(n+1,m))*f(m)+imag(Yb(n+1,m)*e(m)));
        end
* V4 ]  F$ \& l- k% T        PQX(2*n-1)=Pknow(n+1)-P(n);
        PQX(2*n)=Qknow(n+1)-Q(n);
        I(n)=(P(n)-Q(n)*i)/(e(n+1)-f(n+1)*i);
    end
    for n=1:4
        for m=1:4
            YAKEBI(2*n-1,2*m-1)=-imag(Yb(n+1,m+1))*e(n+1)+real(Yb(n+1,m+1))*f(n+1);%Hnm
- n$ F9 T' H0 t8 {9 E            YAKEBI(2*n-1,2*m)=real(Yb(n+1,n+1))*e(n+1)+imag(Yb(n+1,n+1))*f(n+1);%Nnm
$ A% o0 J6 e5 b7 W7 V; G7 [            YAKEBI(2*n,2*m-1)=-real(Yb(n+1,m+1))*e(n+1)-imag(Yb(n+1,m+1))*f(n+1);%Jnm
6 S+ g6 c: t. ~5 P( I            YAKEBI(2*n,2*m)=-imag(Yb(n+1,m+1))*e(n+1)+real(Yb(n+1,m+1))*f(n+1);%Lnm
: q6 ], L3 r% b6 m( n* e        end
9 s2 Z2 x1 g, h$ M        YAKEBI(2*n-1,2*n-1)=-imag(Yb(n+1,n+1))*e(n+1)+real(Yb(n+1,n+1))*f(n+1)+imag(I(n));%Hnn
- c& K; v+ {9 a, v+ [& B" B        YAKEBI(2*n-1,2*n)=real(Yb(n+1,n+1))*e(n+1)+imag(Yb(n+1,n+1))*f(n+1)+real(I(n));%Nnn
% {; f3 X7 j4 e" L( F8 R        YAKEBI(2*n,2*n-1)=-real(Yb(n+1,n+1))*e(n+1)-imag(Yb(n+1,n+1))*f(n+1)+real(I(n));%Jnn
& ]4 {/ A& ?; Y9 P0 }        YAKEBI(2*n,2*n)=-imag(Yb(n+1,n+1))*e(n+1)+real(Yb(n+1,n+1))*f(n+1)-imag(I(n));%Lnn
    end% build the yakebi matrix
' a2 a) P4 A% \7 C9 i    Dfe=inv(YAKEBI)*PQX; %修正fe
# y4 z; j0 m) V; G    for n=1:4
# Y$ Y' O2 o& j        f(n+1)=f(n+1)+Dfe(2*n-1);
        e(n+1)=e(n+1)+Dfe(2*n);
" \& ^1 p4 y, ?0 a  [  N: H  T  d    end
    f
    e
% H. g: }9 F% v# mend
) d0 C! `: W  o; B
+ \! x2 s3 W! V
运行会出现
" g& V: _: T9 X, z$ K- VWarning: Matrix is singular to working precision.
9 a% u& }! ^9 J; q2 g4 M(Type "warning off MATLAB:singularMatrix" to suppress this warning.)
> In f:\MATLAB6p5\work\yy.m at line 35

我这个程序是按照陈珩那本稳态分析的例子编的，想试试。完全没考虑PV节点
( t7 i. S7 X- E0 o1 m这个程序要怎么改进啊？

wangmuyidong

我没有时间仔细给你检查，但是想告诉你，你的雅克比矩阵的形成太麻烦了，极有可能存在错误，并且你得程序也不通用，雅克比矩阵的形成要仔细分析其各个表达式的关系，可以用矩阵运算直接形成。而没有必要用循环。这样会非常耗时的。