1樓：匿名使用者
function lianxuchaoliu
clear;
clc;
n=9;%節點數；
nl=9;%支路數；
isb=1;%平衡節點號；
pr=0.00001;%誤差精度；
b1=[1 4 0.0576i 0 1.05 1; 4 5 0.
017+0.092i 0.158i 1 0; 5 6 0.
039+0.17i 0.358i 1 0; 3 6 0.
0586i 0 1.05 1; 6 7 0.0119+0.
1008i 0.209i 1 0; 7 8 0.0085+0.
072i 0.149i 1 0; 8 2 0.0625i 0 1.
05 1; 8 9 0.032+0.161i 0.
306i 1 0; 9 4 0.01+0.085i 0.
176i 1 0];
%依次是支路首端；末端，支路阻抗；對地電納；支路變比；折算到哪一側標誌（高壓側為1；低壓側為0）；
b2=[0 0 1.05 1.05 0 1; 1.
63 0 1.05 1.05 0 3; 0.
85 0 1.05 1.05 0 3; 0 0 1 0 0 2; 0 0.
9+0.3i 1 0 0 2; 0 0 1 0 0 2; 0 1+0.35i 1 0 0 2; 0 0 1 0 0 2; 0 1.
25+0.5i 1 0 0 2];
%依次是節點的發電機功率；負荷功率；節點電壓初值；pv節點電壓v給定值；節點無功補償裝置容量；節點分類標號（平衡1；pq2；pv3）；
y=zeros(n);%求導納陣；
for i=1:nl
if b1(i,6)==0
p=b1(i,1);q=b1(i,2);
else
p=b1(i,2);q=b1(i,1);
endy(p,q)=y(p,q)-1./(b1(i,3)*b1(i,5));
y(q,p)=y(p,q);
y(q,q)=y(q,q)+1./(b1(i,3)*b1(i,5)^2)+b1(i,4)./2;
y(p,p)=y(p,p)+1./b1(i,3)+b1(i,4)./2;
end%disp('系統的導納陣為：');
%disp(y);
g=real(y);b=imag(y);
for i=1:n
e(i)=real(b2(i,3));
f(i)=imag(b2(i,3));
v(i)=b2(i,4);
endfor i=1:n
s(i)=b2(i,1)-b2(i,2);
b(i,i)=b(i,i)+b2(i,5);
endp=real(s);q=imag(s);
w=zeros(2*n,1);jac=zeros(2*n);
t=0;
while t==0
for i=1:n
if b2(i,6)~=isb
c=0;d=0;
for j=1:n
c=c+g(i,j)*e(j)-b(i,j)*f(j);
d=d+g(i,j)*f(j)+b(i,j)*e(j);
endif b2(i,6)==2%p，q節點；
w(2*i)=p(i)-e(i)*c-f(i)*d;
w(2*i-1)=q(i)-f(i)*c+e(i)*d;
else if b2(i,6)==3%p，v節點；
w(2*i)=p(i)-e(i)*c-f(i)*d;
w(2*i-1)=v(i)^2-(e(i)^2+f(i)^2);
endendelse
w(2*i-1)=0;
w(2*i)=0;
endend%disp(w);
w1=w(3:2*n);
for i=1:n
for j=1:n
if b2(i,6)~=isb
if b2(i,6)==2%p，q節點；
if j~=i
jac(2*i,2*j-1)=-1*(g(i,j)*e(i)+b(i,j)*f(i));
jac(2*i-1,2*j)=(g(i,j)*e(i)+b(i,j)*f(i));
jac(2*i,2*j)=b(i,j)*e(i)-g(i,j)*f(i);
jac(2*i-1,2*j-1)=b(i,j)*e(i)-g(i,j)*f(i);
else if j==i
m=0;h=0;
for r=1:n
m=m+g(i,r)*e(r)-b(i,r)*f(r);
h=h+g(i,r)*f(r)+b(i,r)*e(r);
endjac(2*i,2*j-1)=-1*m-g(i,i)*e(i)-b(i,i)*f(i);
jac(2*i-1,2*j)=-1*m+g(i,i)*e(i)+b(i,i)*f(i);
jac(2*i,2*j)=-1*h+b(i,i)*e(i)-g(i,i)*f(i);
jac(2*i-1,2*j-1)=h+b(i,i)*e(i)-g(i,i)*f(i);
endendelse if b2(i,6)==3%p，v節點；
if j~=i
jac(2*i,2*j-1)=-1*(g(i,j)*e(i)+b(i,j)*f(i));
jac(2*i-1,2*j)=0;
jac(2*i,2*j)=b(i,j)*e(i)-g(i,j)*f(i);
jac(2*i-1,2*j-1)=0;
else if j==i
m=0;h=0;
for r=1:n
m=m+g(i,r)*e(r)-b(i,r)*f(r);
h=h+g(i,r)*f(r)+b(i,r)*e(r);
endjac(2*i,2*j-1)=-1*m-g(i,i)*e(i)-b(i,i)*f(i);
jac(2*i-1,2*j)=-2*f(i);
jac(2*i,2*j)=-1*h+b(i,i)*e(i)-g(i,i)*f(i);
jac(2*i-1,2*j-1)=-2*e(i);
endendendendelse
jac(2*i-1,2*j-1)=0;
jac(2*i,2*j)=0;
jac(2*i-1,2*j)=0;
jac(2*i,2*j-1)=0;
endendend%disp(jac);
jac1=jac(3:2*n,3:2*n);
for k=1:2*n-2
m=0;
for i=k+1:2*n-2
m=jac1(i,k)./jac1(k,k);
w1(i)=w1(i)-m*w1(k);
for j=k+1:2*n-2
jac1(i,j)=jac1(i,j)-m*jac1(k,j);
endendende(2*n-2)=w1(2*n-2)./jac1(2*n-2,2*n-2);
for i=2*n-3:-1:1
c=0;
for k=i+1:2*n-2
c=c+jac1(i,k)*e(k);
e(i)=(w1(i)-c)./jac1(i,i);
endend%disp(e);
for i=1:n-1
e(i+1)=e(i+1)-e(2*i-1);
f(i+1)=f(i+1)-e(2*i);
end%disp(e);
%disp(f);
for i=1:n-1
b(2*i-1)=abs(e(2*i-1));
b(2*i)=abs(e(2*i));
endkb=max(b);
%disp(kb);
if kb1%以切向量中分量最大（絕對值最大）的狀態變數選定為連續引數；
for i=1:2*n-2
jd(i)=abs(d(i));
endfor i=1:2*n-3
if jd(i)>=jd(i+1)
zd=jd(i);ek=i;
else if jd(i)<=jd(i+1)
zd=jd(i+1);ek=i+1;
endendendfor j=1:2*n-1
jj1(2*n-1,j)=0;
endjj1(2*n-1,ek)=1;
w2=zeros(2*n-1,1);
for i=1:2*n-1
w2(i,1)=0;
endif d(ek)>0
w2(ek,1)=1;
else if d(ek)<0
w2(ek,1)=-1;
endendendfor k=1:2*n-1
m=0;
for i=k+1:2*n-1
m=jj1(i,k)./jj1(k,k);
w2(i)=w2(i)-m*w2(k);
for j=k+1:2*n-1
jj1(i,j)=jj1(i,j)-m*jj1(k,j);
endendendd(2*n-1)=w2(2*n-1)./jj1(2*n-1,2*n-1);
for i=2*n-2:-1:1
c=0;
for k=i+1:2*n-1
c=c+jj1(i,k)*d(k);
d(i)=(w2(i)-c)./jj1(i,i);
endend%disp(d);
for i=1:n-1
e(i+1)=e(i+1)+l*d(2*i-1);
f(i+1)=f(i+1)+l*d(2*i);
endt=t+l*d(2*n-1);
disp(d(2*n-1));
%disp(e);
%disp(f);
%disp(t);
%對預估的近似解進行校正；
tt=0;
while tt==0
for i=1:n
if b2(i,6)~=isb
c=0;d=0;
for j=1:n
c=c+g(i,j)*e(j)-b(i,j)*f(j);
d=d+g(i,j)*f(j)+b(i,j)*e(j);
endif b2(i,6)==2% p,q節點；
if b2(i,2)~=0
wj(2*i)=p(i)-t*real(b2(i,2))-e(i)*c-f(i)*d;
wj(2*i-1)=q(i)-t*imag(b2(i,2))-f(i)*c+e(i)*d;
else
wj(2*i)=p(i)-e(i)*c-f(i)*d;
wj(2*i-1)=q(i)-f(i)*c+e(i)*d;
endelse if b2(i,6)==3%p,v節點；
wj(2*i)=p(i)+t*1.5*1.05-e(i)*c-f(i)*d;
wj(2*i-1)=0;
endendelse
wj(2*i)=0;
wj(2*i-1)=0;
endendwj1=wj(3:2*n);
for i=1:n
for j=1:n
if b2(i,6)~=isb
if b2(i,6)==2%p,q節點；
if j~=i
jjac(2*i,2*j-1)=-1*(g(i,j)*e(i)+b(i,j)*f(i));
jjac(2*i-1,2*j)=(g(i,j)*e(i)+b(i,j)*f(i));
jjac(2*i,2*j)=b(i,j)*e(i)-g(i,j)*f(i);
jjac(2*i-1,2*j-1)=b(i,j)*e(i)-g(i,j)*f(i);
else if j==i
m=0;h=0;
for r=1:n
m=m+g(i,r)*e(r)-b(i,r)*f(r);
h=h+g(i,r)*f(r)+b(i,r)*e(r);
endjjac(2*i,2*j-1)=-1*m-g(i,i)*e(i)-b(i,i)*f(i);
jjac(2*i-1,2*j)=-1*m+g(i,i)*e(i)+b(i,i)*f(i);
jjac(2*i,2*j)=-1*h+b(i,i)*e(i)-g(i,i)*f(i);
jjac(2*i-1,2*j-1)=h+b(i,i)*e(i)-g(i,i)*f(i);
endendelse if b2(i,6)==3%p,v節點；
if j~=i
jjac(2*i,2*j-1)=-1*(g(i,j)*e(i)+b(i,j)*f(i));
jjac(2*i-1,2*j)=0;
jjac(2*i,2*j)=b(i,j)*e(i)-g(i,j)*f(i);
jjac(2*i-1,2*j-1)=0;
else if j==i
m=0;h=0;
for r=1:n
m=m+g(i,r)*e(r)-b(i,r)*f(r);
h=h+g(i,r)*f(r)+b(i,r)*e(r);
endjjac(2*i,2*j-1)=-1*m-g(i,i)*e(i)-b(i,i)*f(i);
jjac(2*i-1,2*j)=-2*f(i);
jjac(2*i,2*j)=-1*h+b(i,i)*e(i)-g(i,i)*f(i);
jjac(2*i-1,2*j-1)=-2*e(i);
endendendendelse
jjac(2*i-1,2*j-1)=0;
jjac(2*i,2*j)=0;
jjac(2*i-1,2*j)=0;
jjac(2*i,2*j-1)=0;
endendendjjac1=jjac(3:2*n,3:2*n);
for k=1:2*n-2
m=0;
for i=k+1:2*n-2
m=jjac1(i,k)./jjac1(k,k);
wj1(i)=wj1(i)-m*wj1(k);
for j=k+1:2*n-2
jjac1(i,j)=jjac1(i,j)-m*jjac1(k,j);
endendende1(2*n-2)=wj1(2*n-2)./jjac1(2*n-2,2*n-2);
for i=2*n-3:-1:1
c=0;
for k=i+1:2*n-2
c=c+jjac1(i,k)*e1(k);
e1(i)=(wj1(i)-c)./jjac1(i,i);
endend%disp(e1);
for i=1:n-1
e(i+1)=e(i+1)-e1(2*i-1);
f(i+1)=f(i+1)-e1(2*i);
endfor i=1:n-1
bx(2*i-1)=abs(e1(2*i-1));
bx(2*i)=abs(e1(2*i));
endkb1=max(bx);
if kb10
ttt=0;ccc=ccc+1;
else
ttt=1;
endend
disp(t);
disp(ccc);
%disp(e);
%disp(f);
for i=1:n
fz1(i)=sqrt(e(i)^2+f(i)^2);
enddisp(fz1);
求matlab程式，給定4組資料x1 2，3，7，6，5 1，5，8，
建立模型y f x1,x2,x3,x4 因為y是關於x1,x2,x3,x4的線性函式，所以有f x1,x2,x3,x4 a0 a1 x1 a2 x2 a3 x3 a4 x4,x1 2，3，7，6，5 x2 1，5，8，6，7 x3 4，7，9，12，10 x4 5，8，10，13，11 y 23，2...
設計matlab程式計算圓域上的二重積分
這個可以用matlab的符號積分或者數值積分解決，下面提供4種方法 1 直角座標系符號積分 syms x y int y int sin pi x 2 y 2 y,sqrt 1 x 2 sqrt 1 x 2 先對y積分 i vpa int int y,x,1,1 i 2.0 2 極座標系符號積分 s...
求簡單的金工實習數控車床簡單程式
funac數控車程式設計如下 o9001 n10 g50 x100 z10 設立座標系，定義對刀點的位置 n20 g00 x16 z2 m03 移到倒角延長線，z 軸2mm 處 n30 g01 u10 w 5 g98 f120 倒3 45 角 n40 z 48 加工 26 外圓 n50 u34 w ...