matlab電力系統(tǒng)快速解耦法潮流計算及短路計算程序

電力系統(tǒng)快速解耦法潮流分析及短路計算

程序設計的基本思想:

（1）由于電力系統(tǒng)潮流分析中要利用到矩陣運算，復數運算，故采用mallab編程。采用文件輸入，將

系統(tǒng)的各個參數以文件的形式輸入，便于程序的通用化。

（2）本程序共有兩個輸入文件，分別為線路參數的文件，和已知的節(jié)點狀態(tài)文件（PQ）

（3）為了使程存不僅僅局限于計算9節(jié)點網絡，在形成節(jié)點導納的函數Yn（）中，利用循環(huán)，找出線路

首節(jié)點中的最大編號數，自動確定節(jié)點導納矩陣的維數。故對于任意n節(jié)點網絡，均可以計算出節(jié)

點導納矩陣

（4）在（3）的前提下，為了使程序支持系統(tǒng)增加節(jié)點，增加負荷等造成的PQ參數改變，或者PQ表的

加長。對程序做了如下優(yōu)化。

首先，程序執(zhí)行的基礎是PQ表中平衡節(jié)點在第一行，接下來是PV節(jié)點，最后是PQ節(jié)點，如果系

統(tǒng)添加節(jié)點，或者刪除節(jié)點，均在PQ表的末端操作，會造成PQ表的順序不足平衡節(jié)點、PV節(jié)點、

PQ節(jié)點的順序。故引入了scqcncing（）函數，其作用就是不論輸入的PQ表是什么順序，在程序

讀入后均按平衡-P"PQ的順序排列。

其次，順序打亂的PQ表必須與支路參數表對應，故在Yn（）函數中加入了兩段循環(huán)體，使之對應

（見相應函數體注釋）

（5）在滿足了上述4個條件后，程序便可以通用化了。當然，由于水平有限，且程序未能由大量數據測

試，故缺陷在所難免，這里僅是做了通用化的嘗試。在本文最末附加了該程序通用化的實例。

二、潮流計算框圖

三.定義相應的函數

1.形成節(jié)點導納矩陣的函數Yn()

functionY=Yn(x,y)

與定義一名為Yn的函數，其功能是自動識別輸入表中節(jié)點的個數，形成相應的節(jié)點導納矩陣

[fid,message]=fopen(x,*r*);與從x文件中讀入支路參數

iffid==-l;%判斷文件是否正確打開

else%為變壓器元件

ifn~=0;

Y(m,m)=Y(m,m)+l/(R(i)+lj*X(i));

Y(m,n)=Y(m,n)-l/(k(i)*(R(i)+lj*X(i)));

Y(n,n)=Y(n,n)+l/(k(i)*k(i)*(R(i)+lj*X(i)));

Y(n,m)=Y(n,m)-l/(k(i)*(R(i)+lj*X(i)));

else

Y(m,m)=Y(m,m)+R(i)+1j*X(i);

end;

end;

end;

maxm=HeadNuniber(1);

先通過下面兩個循環(huán)體，確定輸入表中節(jié)點編號的最大值，及為節(jié)點導納矩陣的維數

fori=l:L;

ifmaxm<=HeadNumber(i);

mcixm=IleadNumber(i);

end;

end;

maxn=EndNumber(1);

fori=l:L;

ifmaxn<=EndNumber(i);

maxn=EndNumber(i);

end;

end;

Y=Y(1:max(maxm,maxn),1:max(maxm,maxn));%形成導納矩陣

2.對不滿足要求的PQ參數表進行排序的函數seqencing()

function[Pointstyle,PointNumbcr,Ps,Qs,Uk,Ok]=seqcncing(y)

與定義名為seqencing的函數，其功能是在系統(tǒng)添加節(jié)點，或輸入的PQ參數的順序不滿足要求時，對PQ參

數表進行重新排序，保證平衡節(jié)點放在第一行，接下來是PV節(jié)點，最后是PQ節(jié)點

[fid,message]=fopen(y,'r);*從y文件中讀入PQ參數

iffid==-l;%判斷文件是否正確打開

error(message);

end;

[Pcintstyle,PointNumber,Ps,Qs,Uk,Ok]=textread(y,?%f%f%f%f%f%f*);

fclose(fid);

L=length(PointNumber);

%通過以下兩個循環(huán)體，完成對PQ輸入表的重新排序，其思想是，在PQ參數之前加入一列Pointstyle用

于標識節(jié)點類型，平衡節(jié)點為0,PV節(jié)點為1,PQ節(jié)點為2,以Pointstyle列為基準進行排序

fori=l:L;

forj=l:L-i;

ifPointstyle(j)>Pointstyle(j+l);

t=Pointstyle(j+l);

Pointstyle(j+1)=Pointstyle(j);

Pointstyle(j)=t;

t=PointNumber(j+1);

PointNumber(j+1)=PointNumber(j);

PointNumber(j)=t;

t=Ps(j+l);

Ps(j+l)=Ps(j);

Ps(j)=t；

t=Qs(j+l);

Qs(j+l)=Qs(j);

Qs(j)=t;

t=Uk(j+l);

Uk(j+l)=Uk(j);

Uk(j)=t;

end;

end;

end;

3、形成解耦算法B'矩陣的函數formBI()

functionBl=formBl(x,y)

與定義名為Bl的函數形成解耦算法中的B'矩陣，得到的B'矩陣用B1表示

[fid,message]=fopen(x,*r*);*從x文件中讀入支路參數

iffid==-l;為判斷文件是否正確打開

error(message);

end;

[HeadPoint,HeadNumber,EndPoint,EndNumber,R,X,B,k]=textread(x,*%s%d%s%d%f%f%f%f*);%將

讀人的參數處理為以列為向量的數組

fclose(fid);

L=length(HeadNumber);

[Pcintstyle,PointNumber,Ps,Qs,Uk,Ok]=seqencing(y);

先調用seqencing函數，引入y文件中的PQ參數

A=FointNumber;

先通過以下兩循環(huán)體，實現PQ參數與支路參數的編號對應

fori=l:L;

forj=l:L;

ifHeadNumber(i)==j;

HeadNumber(i)=A(j);

break;

end;

end;

end;

fori=l:L;

forj=l:L;

ifEndNumber(i)==j;

EndNumber(i)=A(j);

break;

end;

end;

end;

Bl=zeros(L,L);

fori=l:L%以行為單位，通過循環(huán)，用支路參數對Bl進行修改，形成夕矩陣

m=HeadNumber(i);n=EndNumber(i);

Bl(m,m)=Bl(m,m)-l/X(i);

Bl(n,n)=Bl(n,n)-l/X(i);

Bl(m,n)=Bl(m,n)+l/X(i);

Bl(n,m)=Bl(n,m)+1/X(i);

end

maxm=HeadNuinber(1);

fori=l:L;

ifmaxm<=HeadNumber(i);

maxm=HeadNumber(i);

end;

end;

maxn=EndNumber(1);

fori=l:L;

ifmaxn<=EndNumber(i);

maxn=EndNumber(i);

end;

end;

B11(2:max(maxm,maxn),2:max(iiaxm,maxn));%形成B'矩陣

4、形成解耦算法B''矩阱的函數formB11()

functionBll=formBll(x,y)

%定義名為BU的函數形成解耦算法中夕’矩陣，用B11表示從x文件中讀入支路參數確定Y,從y文件中

讀入PQ參數確定B11的維數，即除去平衡節(jié)點和pv節(jié)點，此處要求PQ參數錄入時，將平衡節(jié)點和PQ節(jié)

點放在前排，這一要求在Yn函數中通過seqencing函數己經滿足

Y=Yn(x,y);

B=imag(Y);

[Pcintstyle,PointNumber,Ps,Qs,Uk,Ok]=seqencing(y);

i=l；j=l；

whilePointstyle(i)<=1;

i=i+l;

j=j+l；

end;

Bll=B(j:end,j:end);用形成B''矩陣

5、計算正常情況下系統(tǒng)節(jié)點電壓的函數powerflow()

function[UO,00]=powerflow(x,y)

先定義名為powerflow的函數，利用快速解耦算法來計算正常情況下系統(tǒng)內各個節(jié)點的電壓和相角

[Pointstyle,PointNumber,Ps,Qs,Uk,Ok]=seqencing(y);%調用seqencing函數對PQ參數表進行排序

Y=Yn(x,y);為形成節(jié)點導納矩陣，Yn為n維

Bl=formBl(x,y);先形成解耦算法中的B矩陣，B1為n-l維

Bll=formBU(x,y);%形成解耦算法中的B'矩陣，B'為m維

G=real(Y);%?Y的實部

B=imag(Y);為取Y的虛部

U0二Uk;

00=0k;

L=length(PointNumber);

P二zeros(L,1);

Q=zeros(L,1);

dP=zeros(L,1);

dQ=zeros(L,1);

number=l;

i=l;k=l;

whilePointstyle(i)<=l;與通過k值確定系統(tǒng)中PQ節(jié)點的個數

i=i+l;

k=k+l;

end;

whilenumbcr<100%定義迭代次數上限為100次

fori=2:L;

suml=0;

forj=l:L;

suml=suml+UO(j)*(G(i,j)*cos(00(i)-00(j))+B(i,j)*sin(00(i)-00(j)));%潮流方程,n-l

維

end;

dP(i)=Ps(i)-UO(i)*suml;

end

fori=k:L;

sum2=0;

forj=l:L;

sum2=sum2+U0(j)*(G(i,j)*sin(00(i)-00(j))-B(i,j)*cos(00(i)-00(j)));%潮流方程，m

維

end;

dQ(i)=Qs(i)-U0(i)*sum2;

end

dPl=dP(2:L)./U0(2:L);

dQl=dQ(k:L)./UO(k:L);

a=max(norm(dPl,inf));

b-max(norm(dQl,inf));

ifmax(a,b)<0.00001%判斷是否收斂

break;

disp('迭代’)

disp(k);

disp('次后收斂');

else%如不收斂，

dO=-inv(Bl)*dPl;%dO為nT維

dU=-inv(Bll)*dQl;斷iU為m維

zerol=zeros(k-1,1);

zcro2=[0];

DU=[zerol;dU];

DO=[zero2;dO];

UO=UO+DU;

00=00+DO;

number=number+1;

end;

ifnumber==100;

dispC迭代100次后不收斂，迭代結束')；

end;

end;

四.對相應的系統(tǒng)進行潮流分析和短路計算

定義完上述函數之后，可宜接調用函數形成導納矩陣，計算正常情況下的節(jié)點電壓，進行短路計算計算短

路也流，短路后各個節(jié)點也壓以及支路潮流分布。

程序的輸入表共有三個,為4network,txt*,'network2.txt','PQ參數.txt'。其中'network.txt'

存放不含發(fā)電機節(jié)點導納和負荷導納的系統(tǒng)支路參數，'network2.txt'存放含發(fā)電機節(jié)點導納和負荷導納

的系統(tǒng)支路參數，'PQ參數.txt'存放系統(tǒng)的PQ參數，現列表如下

network-記事本

文件(F)期(E)颯O)X(v)幫助(H)

BUS1BUS400.05761.01.0

BUS2BUS700.06251.01.0

BUS3BUS900.05861.01.0

BUS4BUS50.010.0850.0880

BUS4BUS60.0170.0920.0790

BUS5BUS70.0320.1610.1530

BUS6BUS90.0390.170.1790

BUS7BUS80.00850.0720.07450

BUS8BUS90.01190.10080.10450

第一列:首節(jié)點名稱第二列:首節(jié)點編號第三列:尾節(jié)點名稱第四歹U：尾節(jié)點編號

第五列:支路參數R第六列:支路參數X第七列:支路參數B/2第八列:變比K

\PQ參數-記事本

文件(F)嗝號(E)悟式(O)直言英超麗彳

01

12001.040

131.6301.0250

240.8501.0250

250010

26-1.25-0.510

27-0.9-0.310

280010

29-1-0.3510

0010

第一列:節(jié)點類型,o-平衡節(jié)點,i-py節(jié)點,2-PQ節(jié)點第二列:節(jié)點編號

第三列:節(jié)點的有功P第四列:節(jié)點的無功Q第五列:節(jié)點電壓,PQ節(jié)點置為0

第六列:節(jié)點電壓的相角,初始值為0

netv/ork2-記事本

文件(F)編輯(E)格式(O)查看(V)幫助(H)

11

BUSIBUS4q0576.01.0

21.

BUS23BUS770625101.0

c.

BUS34BUS9y0586O.01.0

BUS44BUS550.01085O.0880

BUS4BUS660.017092.0790

5O.

BUSS6BUS?70.032161O.1530

BUS67BUS9c0.03917O.1790

.

BUS78BUS8g0.0085072O07450

o

BUS8BUS9c0.0119100810450

1Oo

BUSI0-3.33333

2OcO

BUS23O0c-3.33333

BUS34O0c-3.33333o

BUS45O0c0o

BUSS6O0-1.26110.5044o

BUS67O0-0.87760.2925o

BUS78O0C0o

BUS89O0-0.96890.3391O

BUS90C0

程序代碼如下:

%正常情況下節(jié)點導納矩陣

Y=Yn('network,txt，，PQ參數.txt');

disp('節(jié)點導納矩陣')

disp(sparse(Y));

為正常情況下潮流計算結果兵繪制潮流計算框圖

[U0,00]=powerflow(*network,txt*,*PQ參數.txt')%各個節(jié)點電壓為U0,相角為00

disp('正常情況下系統(tǒng)各個節(jié)點電壓');

disp(UO);

disp('正常情況下系統(tǒng)各個節(jié)點相角，);

disp(00);

00=00*180/pi

disp('正常情況下系統(tǒng)各個節(jié)點相角角度值’)

%計算平衡節(jié)點功率和PV節(jié)點無功功率：

form=l:9

pl=pl+U0(1)*U0(m)*(G(1,m)*cos(00(1)-00(m))+B(l,m)*sin(00(1)-00(m)));

ql=ql+U0(l)*U0(in)*(G(1,m)*sin(00(l)-OO(ni))-B(l,m)*cos(00(1)-00(m)));

q2=q2+U0(2)*U0(m)*(G(2,m)*sin(00(2)-00(m))-B(2,m)*cos(00(2)-00(m)));

q3=q3+U0(3)*U0(m)*(G(3,m)*sin(00(3)-00(m))-B(3,m)*cos(00(3)-00(m)));

end;

disp('平衡節(jié)點有功、無功，兩個電壓節(jié)點有功分別為')

disp(pl);%有功功率和無功功率

disp(ql);%平衡節(jié)點無功功率

disp(q2);%PV節(jié)點2無功功率

disp(q3);%PV節(jié)點3無功功率

為增加發(fā)電機導納yi和負荷導納yLDi后的發(fā)電機節(jié)點和負荷節(jié)點的自導納

Y=YnCnetwork2.txtJ,*PQ參數.txt')%network2文件口存放的是經修改的支路參數

disp('加發(fā)電機導納yi和負荷導納yLDi后的發(fā)電機節(jié)點和負荷節(jié)點的自導納');

diag(Y)

用精確法計算短路電流，短路后各節(jié)點電壓以及網絡中各支路的電流分布

Y=Yn(，network?.Ixl','PQ參數.Ixl');先形成包括發(fā)電機內阻抗和符合阻抗節(jié)點導納矩陣

Z=inv(Y);%形成節(jié)點阻抗矩陣％題中給出的是4節(jié)點短路

Zf=Z(:,4)；

[U0,00]=powerflow(，network,txt','PQ參數.txt');

U01=U0.*cos(00)+lj*UO.*sin(OO);

If=U01(4)/Z(4,4);%計算短路電流If

U=U01-Zf*If;斷十算短路后的各個節(jié)點電壓

[fid,message]=fopen(，network2.txt*,*r*);*從x文件中讀入支路參數

iffid==-l;%判斷文件是否正確打開

error(message);

end;

[UeadPoint,HeadNumber,EndPoint,EndNumber,X,B,k]=textreadnetwork2.txt','%s%s%d%f

%f%f%f');

row=[HoadNumber,EndNumber];

fori=l:9%計算系統(tǒng)中各個支路的潮流分布

I(row(i,1),row(i,2))=(U(row(i,1))-U(row(i,2)))/Z(row(i,1),row(i,2));

end;

disp('阻抗矩陣中的第f列');

disp('精確算法短路電流If');

disp('If模值’);

disp(abs(If));

disp('相角為');

disp(angle(If));

disp('短路后各個節(jié)點的電壓');

disp(U);

disp('模值為');

disp('短路后各個支路的電流值I');

disp('模值為r);

disp('相角為(度)');

與近似計算短路電流，短路后各節(jié)點電壓以及網絡中各支路的電流分布

If2=l/Z(4,4);%計算短路電流If

U2=l-Zf*If;先計算短路后的各個節(jié)點電壓

disp('近似算法短路電流U2f);

disp('If2模值');

disp(abs(If));

disp('相角為(度)')；

disp(angle(If2));

disp('短路后各個節(jié)點的電壓');

disp(U2);

disp('模值為');

程序運行結果如下:

節(jié)點導納矩陣

Y=

(1,1)0-17.36111

(4,1)0+17.361li

⑵2)0-16.OOOOi

(7,2)0+16.OOOOi

(3,3)0-17.0648i

⑼3)0+17.0648i

(1,4)0il7.3611i

(4,4)3.3074-39.3089i

(5,4)-1.3652+11.6041i

(6,4)-1.9422+10.5107i

(4,5)-1.3652+11.6041i

⑸5)2.5528-17.3382i

(7,5)-1.1876+5.9751i

(4,6)-1.9422+10.5107i

(6,6)3.2242-15.8409i

(9,6)-1.2820+5.58821

(2,7)0+16.OOOOi

⑸7)-1.1876+5.9751i

(7,7)2.8047-35.4456i

(8,7)-1.6171+13.6980i

(7,8)-1.6171+13.6980i

(8,8)2.7722-23.3032i

⑼8)-1.1551+9.7843i

(3,9)0+17.0648i

(6,9)-1.2820+5.5882i

(8,9)-1.1551+9.7843i

⑼9)2.4371-32.1539i

迭代k=9

次后收斂

正常情況下系統(tǒng)各個節(jié)點電E

U0=

1.0400

1.0250

1.0250

1.0258

0.9956

1.0127

1.0258

1.0159

1.0324

正常情況下系統(tǒng)各個節(jié)點相角

00=

0

0.1620

0.0814

-0.0387

-0.0696

-0.0644

0.0649

0.0127

0.0343

正常情況下系統(tǒng)各個節(jié)點相角角度值

00=

0

9.2800

4.6647

-2.2168

-3.9888

-3.6875

3.7197

0.7275

1.9667

平衡節(jié)點有功、無功，兩個電壓節(jié)點有功分別為

pl=

0.7164

ql=

0.2705

q2=

0.0665

q3=

-0.1086

加發(fā)電機導納yi和負荷導納yLDi后的發(fā)電機節(jié)點和負荷節(jié)點的自導納

ans=

0-20.6944i

0-19.3333i

0-20.3982i

3.3074-39.3089i

1.2917-16.8338i

2.3466-15.5484i

2.8047-35.4456i

1.8033-22.964li

2.4371-32.1539i

阻抗矩陣中的第f列

Zf=

-0.0630+0.1475i

-0.0645+0.0860i

-0.0634+0.0864i

-0.0751+0.1758i

-0.0892+0.1534i

-0.0854+0.1532i

-0.0779+0.1040i

-0.0832+0.1032i

-0.0758+0.1032i

精確算法短路電流

If=

-2.2963-4.8486i

If模值為

5.3649

相角為

-115.3422

短路后各個節(jié)點的電壓

U=

0.1801+0.0333i

0.4463+0.0502i

0.45730.0257i

-0.00000.0000i

0.04460.1498i

0.07160.1271i

0.3406-0.0725i

0.3243-0.1533i

0.3572-0.0950i

模值為

0.1831

0.4491

0.4580

0.0000

0.1563

0.1459

0.3482

0.3587

0.3696

短路后各個支路的電流值I

(1,4)0.25011.1140i

(4,5)0.8560—0.2072i

(4,6)0.8320-1-0.004li

(2,7)0.5086—0.9694i

(5,7)0.681811.8688i

(7,8)0.3965—0.3175i

(3,9)0.2141—0.7887i

(6,9)0.8241+1.7182i

(8,9)-0.235510.3676i

模值為I

(1,4)1.1417

(4,5)0.8807

(4,6)0.8320

(2,7)1.0947

(5,7)1.9893

(7,8)0.5079

(3,9)0.8172

(6,9)1.9056

(8,9)0.4365

相角為（度）

(1,4)-102.6530

(4,5)-13.6077

(4,6)0.2855

(2,7)-62.3148

(5,7)69.9555

(7,8)-38.6844

(3,9)-74.8141

(6,9)64.3752

(8,9)122.6426

近似算法短路電流

If2=

0.2105-5.2357i

If2模值5.2399

相角為(度)-87.6981

短路后各個節(jié)點電壓值

U2=

0.1611?0.OOOOi

0.4785-0.0351i

0.4773-0.03841

0.0000

0.1155-0.0063i

0.1178-0.0098i

0.3698-0.0424i

0.3673-0.0436i

0.3752-0.0459i

模值為

0.1611

0.4797

0.4788

0.0000

0.1157

0.1182

0.3722

0.3698

0.3780

精確算法與近似算法的比較分析

（1）短路電流誤差

（If精Tf近）/If精=2.33%

\4）精-4>近\=27.6441（度）

兩電流相差進2.33%,可見在精度要求不是很高的前提下，近似算法是完全可取的。

（2）電壓模值誤差比較

（U精-U近）川精二

0.1205

-0.0681

-0.0455

0.0000

0.2597

0.1900

-0.0690

-0.0311

-0.0227

可見電壓模值的誤差均在25與以內。

至此，上述程序已完成題目的所有要求。為了實現程序的通用行，在編寫程序的時候對程序進行了優(yōu)化，

使得上述程序不但能夠滿足本題要求，而且對于任意多節(jié)點的系統(tǒng)，只要將支路參數，PQ參數列表輸入ixl

文件中，有上述程序均可完成正常情況下潮流計算和短路計算。

而且，上述程序還支持支路參數的修改，添加或刪除支路、節(jié)點。

四.附加：程序通用化實例

若在原題的基礎上，在8母線和9母線之間增加一條母線10,,相應的支路參數如下表

支路R(p.u.)X(p.u.)B/2(L)(p.u.)

「400.05761.0

2~700.06251.0

3~900.05861.0

4~50.010.G850.088

4飛0.0170.G920.079

5~70.0320.1610.153

6~90.0390.170.179

7~80.00850.0720.0745

8?100.00590.05080.1045

9~100.00500.00500.1045

新增10節(jié)點為PQ節(jié)點P=0,Q=0；

由上述條件可以寫出添加節(jié)點后的支路參數文件為

]1-記事本

文件(F)編號(E)嗓(O)M(V)幫助(H)

US114o

US22BUS47o0.05761.01.0

US33BUST9o0.06251.01.0

US44BUS950.05861.01.0

-n

US44BUS56010.0850.0880

US55BUS670170.0920.0790

US66BUST90320.1610.1530

US77BUS980390.170.1790

US88BUS81O00850.0720.07450

US89BUS91300590.05080.10450

BUS900600.05000.10450

PQ參數文件為

32-記事本

文件(F)編輯(E)格式(O)查看(V)幫助而

o

0101.04O

126301.025O

138501.025

24010o

25.25-0.51O

26.9-0.31

27010O

28-0.35o

2901O

21001

以這兩個表作為輸入，執(zhí)行上述程序，計算潮流，得到各個節(jié)點的電壓和相角為

U0=

1.0400

1.0250

1.0250

1.0278

0.9988

1.0168

1.0304

1.0265

1.0397

1.0385

00=

0

0.1601

0.0795

-0.0386

-0.0695

-0.0643

0.0635

0.0113

0.0327

0.0215

即可得到十個節(jié)點的電壓和相角，完成了程序的通用化。

五.結語

學習是一個理論與實踐相結合的過程，學習電力系統(tǒng)分析，對相應的電力系統(tǒng)進行簡單的潮流計算和短路

分析是我們必須掌握的基本知識。在用計算機編程是手段實現系統(tǒng)分析的功能時，我們的問題是只著眼于

既定的題目和系統(tǒng)。所以，怎樣將程序通用化是?一個值得考慮的問題，由于電力系統(tǒng)的復雜性以及個人知

識水平的限制，謹在此方面做了一-些改進，真正的知識體系還亟待在今后的學習中得以建立和提高。附錄

資料：matlab畫二次曲面

一、螺旋線

1.靜態(tài)螺旋線

a=0:0.1:20*pi;

,,

h=plot3(a.*cos(a),a.*sin(a)z2.*a,b/linewidth\2);

axis([-50,50,-50,50,0,150]);

gridon

set：h,'erasemode','none','markersize',22);

xlabel('x軸');ylabel('y^');zlabel('z釉')；

title('靜態(tài)螺旋線');

2.動態(tài)螺旋線

t=0:0.1:10*pi;

i=l；

,,,

h=plot3(sin(t(i)),cos(t(i)),t(i)/*',erasemode'z'none);

gridon

axis([-22-22035])

fori=2:length(t)

,,,

set(h,xdata/sin(t(i))/ydata'/cos(t(i))/zdata'/t(i));

drawnow

pause(O.Ol)

end

titled動態(tài)螺旋線'）;

（圖略）

3.圓柱螺旋線

t=0:0.1:10*pi;

x=r.*cos(t);

y=r.*sin(t);

z=t;

plot3(xzy,z/h','linewidth',2);

gridon

axis('square')

xlabel('x軸)ylabel('y軸');zlabel('z軸');

title。圓柱螺旋線')

圖樣守亶毆

二、旋轉拋物面

b=0:0.2:2*pi;

rx,Yl=meshqrid(-6:0.1:6);

Z=[X.人2+Y.人2)./4;

meshc(X,Y,Z);

axis('square')

xlabel('x軸)ylabel('y軸');zlabel('z軸)

shadingflat;

title('旋轉拋物面')

或直接用：ezsurfc('(X.A2+Y.A2)./4')

三、橢圓柱面

loadclown

,,

ezsurf('(2*cos(u))'/4*sin(u)/'v'/[0,2*pi/0z2*pi])

%視角處理

view(-105z40)

shadinginterp%燈光處理

colormap(map)%顏色處理

gridon%添加網格線

axisequal%使x,y軸比例一致

xlabel('x軸');ylabel('y軸');zlabel('z軸');

shadingflat;

title。橢圓柱面')％添加標題

松回住囹

四、橢圓拋物面

b=0:0,2:2*pi;

[X,Y]=meshgrid(-6:0.1:6);

Z=X.人2./9+Y.人2./4;

meshc(X,Y,Z);

axis('square')

xlabel('x軸');ylabel('y軸)zlabel(2軸');

shadingflat;

title('橢圓拋物面’)

或直接用：ezsurfc('X.A2./9+Y.A2./4')

楙世糖州直

五、'雙葉雙曲面

,

ezsurf('8*tan(u)*cos(v)/'8.*tan(u)*sin(v)'z'2.*sec(u)',[-pi./2/3*pi./2,0,2*pi])

axisequal

gridon

axissquare

xlabel('x?i');ylabel('y軸');zlabel('z軸');

shadingflat;

title（'雙葉雙曲面'）

六、雙曲柱面

loadclown

,,,,

ezsurf('2*sec(u)'/2*tan(u)zvz[-pi/2/pi/2,-3*piz3*pi])

holdon%在原來的圖上繼續(xù)作圖

ezsurf('2*sec(u)72*tan(u)7v'z[pi/2/3*pi/2z-3*piz3*pi])

colormap(map)

shadinginterp

view(-15z30)

axisequal

gridon

axisequal

xlabel('x軸');ylabel('y軸');zlabel('z軸');

shadingflat;

title('雙曲柱面')

僅6線面

七、雙曲拋物面(馬鞍面)

[X,Y]=meshgrid(-7:0.1:7);

Z=X.A2./8-Y.A2./6；

meshc(X,Y,Z);

view(85,20)

axis('square')

xlabel('x軸');ylabel('y