第一篇：潮流計算--電力系統大作業(C++編寫)
程序設計所涉及二叉樹解釋
任意一棵樹或一個森林都能唯一地對應一棵二叉樹，由此而編寫本程序。本程序采用類二叉樹為整體結構，二叉樹類下定義節點類，每一條支路均為樹的一個節點，支路所有的參數均作為節點的屬性，并給節點加入屬性“支路編號”，并以支路編號為依據構建二叉樹，這就要求提前根據二叉樹結構給每一個支路編號。
支路編號原則：左子樹上所有編號均小于其雙親的編號，右子樹上所有編號均大于其雙親的編號，為了便于查看，本程序在節點較少時編號從1開始，逐個遞加至支路數；當支路較多時，可不必拘泥于逐個遞加，只要滿足支路編號原則即可。
例如習題3-4：
123118kV113kV8.5+j20.5Wj2.82×10-4S1.22+j20.2Wj2.82×10-4SDST=0.17+j1.7MVA40+j30MVA20+j15MVA
程序二叉樹結構示意圖：
本二叉樹中節點1即支路1為題目中節點1與節點2之間的部分； 本二叉樹中節點2即支路2為題目中節點2與節點3之間的部分； 對于習題3-3：
程序二叉樹結構示意圖：
本二叉樹中節點2即支路2為題目中節點1與節點2之間的部分； 本二叉樹中節點1即支路1為題目中節點2與節點3之間的部分； 本二叉樹中節點3即支路3為題目中節點2與節點4之間的部分.拓展：
如下多支路網絡：
對于三節點網絡需先進行以下網絡處理轉化為標準二叉樹，而后進行計算；
圖中三角形表明該項阻抗為零，為純導線，并進行相應參數補充進行計算。
程序說明文檔
******************************************************************************************************************************************************* 本程序測試使用方法：在E盤根目錄下建立輸入文件：
輸入文件名：input.txt;
將所附算例對應輸入文件內容復制粘貼至上述文件中，在VC++6.0環境下運行cpp文件得出結果。
輸出文件在E盤根目錄下
輸出文件名為：data.txt;******************************************************************************************************************************************************* 程序功能說明：本程序可以計算任意長度線型開始網絡潮流；
支持多電壓等級下的計算；
可在除供電節點外任意節點引出負載；
*******************************************************************************************************************************************************
輸入格式說明：以支路為基本單位，按潮流方向輸入數據：
以下例示意：
//支路個數
1,110,118,8.5,20.5,0.000564,0,0,1,0,0
//線路等效的支路
2,11,110,1.22,20.2,0,40,30,10,0.17,1.7
//變壓器等效的支路
.//按此方式知道輸入所有的支路
.//輸入從上到下的順序為潮流在
.//線型開式網絡中的流動方向
/*
第一行輸入支路個數，回車
第二行至后輸入各個支路參數，回車分隔不同支路；
各行輸入的支路參數順序是：
支路編號，末端電壓，始端電壓，線路等效電阻，線路等效感抗，線路等效容納，末端輸入有功，末端輸入無功，變比，變壓器有功勵磁損耗，變壓器無功勵磁損耗
對于線路等效電路：變壓器有功/無功損耗輸入零，變比輸入1；
對于變壓器等效電路：所有的參數均歸算至高壓側，Rt,Xt對應輸入線路等效電阻/感抗的位置，線路等效容抗為零；
對于個節點的引出負荷：輸入至以此節點為末節點的支路的末端輸入有功/無功部分；
*/ ******************************************************************************************************************************************************** 輸出文件格式說明：輸出的內容包括
（1）支路信息：每個支路元件的始端有功、無功和末端有功、無功；有功損耗無功損耗；電壓損耗；
（2）全網信息：全網的總電源有功、總負荷有功、有功損耗、網損率；
（3）迭代信息：每次完整迭代后的所有內容；
具體在輸出文件中都明確標出。
********************************************************************************************************************************************************* 程序中變量定義說明：
類中定義的變量
class line_part{
//定義支路類
double U[2];
//支路電壓降落：U[0]電壓降落橫分量，U[1]電壓降落縱分量
double k;
//變壓器變比
double val;
//支路排序
double U_end;
//支路末端電壓
double U_begin;
//支路首段電壓
double X[3];
//支路等效阻抗: X[0]電阻，X[1]感抗，X[2]容納
double S_end[2];
//支路末端功率：S_end[0]有功，S_end[1]無功
double S_begin[2];
//支路首段功率：S_begin[0]有功，S_begin[1]無功
double S0[2];
//變壓器勵磁損耗：S0[0]有功，S0[1]無功
double S_org[2];
//支路末端負載：S_org[0]有功，S_org[1]無功
line_part *lchild,*rchild;
//支路的后繼兩個支路 }
class BinTree{
//定義樹類
void PreOrder(){PreOrder(root);};
//樹的先序遍歷修改電壓
void PostOrder(){PostOrder(root);};
//樹的后序遍歷修改潮流
void display(){display(root);};
//樹的先序遍歷顯示數據
line_part *root;
//樹的根，是一個支路類
};主函數中定義的數據
ofstream outfile;
//輸出數據流定
ifstream infile;
//輸入數據流定
const int M(a);
//支路個數常量
主程序（復制粘貼到C++就能用）
**************************************************************************************************************************************************************
#include
double p_cost_all=0;
double sqr(double x){
//平方計算函數 return x*x;};
class line_part{
//定義支路類 private:
double val;
//支路排序
double U_end;
//支路末端電壓
double U_begin;
//支路首段電壓
double X[3];
//支路等效阻抗: X[0]電阻，X[1]感抗，X[2]容納
double S_end[2];
//支路末端功率：S_end[0]有功，S_end[1]無功
double S_begin[2];
//支路首段功率：S_begin[0]有功，S_begin[1]無功
double S0[2];
//變壓器勵磁損耗：S0[0]有功，S0[1]無功
double S_org[2];
//支路末端負載：S_org[0]有功，S_org[1]無功
line_part *lchild,*rchild;public:
double U[2];
//支路電壓降落：U[0]電壓降落橫分量，U[1]電壓降落縱分量
double k;
//變壓器變比
public: line_part(){ val=0;U_end=0;U_begin=0;X[0]=0;X[1]=0;
X[2]=0;S_end[0]=S_org[0]=0;S_end[1]=S_org[1]=0;
S_begin[0]=0;S_begin[1]=0;k=1;
S0[0]=0;S0[1]=0;U[0]=0;U[1]=0;
lchild = rchild= NULL;}
line_part(double vall,double u_end=0,double u_begin=0,double r=0，double x=0,double b=0,double Pe=0,double Xe=0,double K=0,double P0=0,double Q0=0){
val=vall;
U_end=u_end;
U_begin=u_begin;
X[0]=r;
X[1]=x;
X[2]=b;
S_end[0]=S_org[0]=Pe;
S_end[1]=S_org[1]=Xe;
S_begin[0]=0;
S_begin[1]=0;
k=K;
S0[0]=P0;
S0[1]=Q0;
U[0]=0;
U[1]=0;
lchild = rchild= NULL;
};friend class BinTree;
friend void pass_U(line_part*a,line_part*b,line_part*c);
//電壓傳遞函數
friend void pass_w(line_part*a,line_part*b,line_part*c);
//功率傳遞函數
friend void pass_U2(line_part*a,line_part*b);
//電壓傳遞函數
friend void pass_w2(line_part*a,line_part*b);
//功率傳遞函數
void Sbegin(){
(U_end)*X[2]/2))/sqr(Uend);;-sqr(U_begin)*X[2]/2-sqr(U_end)*X[2]/2;};void Uend(){
double U_heng(0),U_zong(0);
X[2]/2-S0[1];1])/U_begin;
//U_heng即是△u2 0])/U_begin;
//U_zong即是δu2 sqr(U_zong))/k;
double get_val(){
};
//支路首段功率計算函數 double Uend,I2;Uend=k*U_end;I2=(sqr(S_end[0])+sqr(S_end[1]-sqrS_begin[0]=S_end[0]+I2*X[0]+S0[0]S_begin[1]=S_end[1]+I2*X[1]+S0[1]
//支路末端電壓計算函數 double p_begin,q_begin;p_begin=S_begin[0]-S0[0];q_begin=S_begin[1]+sqr(U_begin)*U_heng=(p_begin*X[0]+q_begin*X[U_zong=(p_begin*X[1]-q_begin*X[U_end=sqrt(sqr(U_begin-U_heng)+U[0]=U_heng;U[1]=U_zong;};
//返回支路編號 if(this==0){return-1;}else{ if(val>0&&val<100){return val;} else return-1;}
double get_Uend(){
//返回支路末端電壓
return U_end;};
double get_Ubegin(){
//返回支路首段電壓
return U_begin;};
double get_Pbegin(){
//返回支路首段有功
return S_begin[0];};
double get_Pend(){
//返回支路末端有功
return S_end[0];};
double get_Qbegin(){
//返回支路首段無功
return S_begin[1];};
double get_Qend(){
//返回支路末端無功
return S_end[1];};
double get_Pcost(){
//返回支路有功損耗
return S_begin[0]-S_end[0];};
double get_Qcost(){
//返回支路無功損耗
return S_begin[1]-S_end[1];};line_part *get_lchild(){
//返回支路無功損耗
return lchild;};line_part *get_rchild(){
//返回支路無功損耗
return rchild;};};
void pass_U(line_part *a,line_part *b,line_part *c){(*c).U_begin=(*a).U_end;(*b).U_begin=(*a).U_end;};void pass_w(line_part *a,line_part *b,line_part *c){(*a).S_end[0]=(*b).S_begin[0]+(*c).S_begin[0]+(*a).S_org[0];(*a).S_end[1]=(*b).S_begin[1]+(*c).S_begin[1]+(*a).S_org[1];};void pass_U2(line_part *a,line_part *b){(*b).U_begin=(*a).U_end;};void pass_w2(line_part *a,line_part *b){(*a).S_end[0]=(*b).S_begin[0]+(*a).S_org[0];(*a).S_end[1]=(*b).S_begin[1]+(*a).S_org[1];};
class BinTree{ public:
friend void pass_U(line_part*,line_part*,line_part*);
//電壓傳遞函數
pass_w(line_part*,line_part*,line_part*);
line_part(0,0,0,0,0,0,0,0,0,0,0);root = aa;}
vall,double u_end=0,double u_begin=0,double r=0, Pe=0,double Xe=0, P0=0,double Q0=0){ egin,r,x,b,Pe,Xe,K,P0,Q0);
private:
*&t,double vall,double u_end,double u_begin,double r, Pe,double Xe, Q0);
};
friend void //功率傳遞函數 BinTree(){line_part *aa=new line_part *Getroot(){return root;} void insertline_part(double
double x=0,double b=0,double
double K=0,double insertline_part(root,vall,u_end,u_b} void PreOrder(){PreOrder(root);};void PostOrder(){PostOrder(root);};void display(){display(root);};line_part *root;void insertline_part(line_part
double x,double b,double
double K,double P0,double void PreOrder(line_part *&t);void PostOrder(line_part *t);void display(line_part *&t);
void BinTree::insertline_part(line_part *&t, double vall,double u_end=0,double u_begin=0,double r=0，double x=0,double b=0,double Pe=0,double Xe=0，double K=0,double P0=0,double Q0=0){
//插入節點
double ass=t->get_val();if(t==0||t->get_val()<=0){
t=new line_part(vall,u_end,u_begin,r,x,b,Pe,Xe,K,P0,Q0);} else if(vall
insertline_part(t->lchild, vall,u_end,u_begin,r,x,b,Pe,Xe,K,P0,Q0);} else {
insertline_part(t->rchild, vall,u_end,u_begin,r,x,b,Pe,Xe,K,P0,Q0);} };/////////////////////////////////////////////////////////////////////////////////////////////////////////////// void BinTree::display(line_part *&t){
if(t->get_val()>0&&t->get_val()<100){
display(t->lchild);
display(t->rchild);
ofstream outfile1;
outfile1.open(“e:data.txt”,ios::ate);
p_cost_all +=t->get_Pcost();
double U;
// 計算并存放各個節點的電壓相角（始端為零）
U=atan2(t->U[1],(t->get_Ubegin()-t->U[0]))/3.1415926*180;
outfile1<<“支路”<
”< outfile1<<“ ”<<“首端電壓”<<“ ”< outfile1<<“ ”<<“末端電壓”<<“ ”< //計算并存放各個支double U_cost; //計算并存放各個支路的電壓損耗 U_cost=t->get_Ubegin()-t->get_Uend(); outfile1<<“ ”<<“電壓損耗”<<“ outfile1<<” “<<”始端有功“<<” 首端電壓 outfile1<<“ ”<<“始端無功”<<“ 路的末端電壓 outfile1<<” “<<”末端有功“<<” 首端電壓 outfile1<<“ ”<<“末端無功”<<“ 路的末端電壓 double P_cost; 的有功損耗 P_cost=t->get_Pcost(); outfile1<<” “<<”有功損耗“<<” double Q_cost; 的無功損耗 Q_cost=t->get_Qcost(); outfile1<<“ ”<<“無功損耗”<<“ outfile1.close(); ”< //計算并存放各個支“< //計算并存放各個支 //計算并存放各個支路 “< //計算并存放各個支路 ”< };///////////////////////////////////////////////////////////////////////////////////////////////////////////////// void BinTree::PreOrder(line_part *&t){ //先序遍歷二叉樹修改電壓 0) get_lchild()->get_val()<100) get_rchild()->get_val()<100) get_rchild()->get_val()<100) if(t->get_val()>0&&t->get_val()<10{ if(t->get_lchild()->get_val()>0&&t-> { if(t->get_rchild()->get_val()>0&&t-> { t->Uend(); pass_U(t,t->lchild,t->rchild); t->lchild->Uend(); t->rchild->Uend(); } else{ t->Uend(); pass_U2(t,t->lchild); t->lchild->Uend();} } else { if(t->get_rchild()->get_val()>0&&t-> { t->Uend(); pass_U2(t,t->rchild); t->rchild->Uend(); } else{} } PreOrder(t->lchild);PreOrder(t->rchild); } };void BinTree::PostOrder(line_part *t){ //后序遍歷二叉樹修改潮流 if(t->get_val()>0&&t->get_val()<100) get_lchild()->get_val()<100) get_rchild()->get_val()<100) get_rchild()->get_val()<100) { PostOrder(t->lchild);PostOrder(t->rchild);if(t->get_lchild()->get_val()>0&&t-> { if(t->get_rchild()->get_val()>0&&t-> { t->lchild->Sbegin(); t->rchild->Sbegin(); pass_w(t,t->lchild,t->rchild); t->Sbegin(); } else{ t->lchild->Sbegin(); pass_w2(t,t->lchild); t->Sbegin();} } else { if(t->get_rchild()->get_val()>0&&t-> { t->rchild->Sbegin(); pass_w2(t,t->rchild); t->Sbegin(); } else{} }; void main(){ ofstream outfile; //輸入數據流定義 infile>>a; cout<<“節點個數 ”< const int M(a); 量 double *A=new double[10*M];支路數據數組 while(l<11*M){ } } //輸出數據流定義 outfile.open(“e:data.txt”);outfile.clear();outfile.close();ifstream infile;infile.open(“e:input.txt”); int a;int l(0);char b; //支路節點數常 //輸入流輸入各 infile>>A[l]; infile.get(b);l++;};outfile.open(“e:data.txt”,ios::ate);outfile<<“節點個數 ”< i=0;i *i+1],A[11*i+2],A[11*i+3],A[11*i+4]，*i+8],A[11*i+9],A[11*i+10]); 如下“< outfile.close(); BinTree elec;for(int //添加節點 elec.insertline_part(A[11*i+0],A[11 A[11*i+5],A[11*i+6],A[11*i+7],A[11} for(i = 0;i<5;i++){ outfile.open(”e:data.txt“,ios::ate);outfile< outfile< elec.PostOrder(); elec.PreOrder();elec.display(); outfile.open(“e:data.txt”,ios::ate); “< ”< “< ”< outfile< line_part *ee=new line_part();ee=elec.Getroot(); double aa =ee->get_Pbegin(); outfile<<“全網的總電源有功 outfile.close(); outfile.open(”e:data.txt“,ios::ate);double bb;bb =aa-p_cost_all; outfile<<”全網的總負荷有功 outfile<<“全網的總有功損耗 double cc; cc = p_cost_all/aa*100; outfile<<”全網的網損率 outfile.close(); p_cost_all=0.0; } }; **************************************************************************************************************************************************************習題3-3 input: 3 2,35,38.5,1.2,2.4,0,3,2,1,0,0 1,35,35,1,2,0,5,3,1,0,0 3,35,35,2,4,0,2,3,1,0,0 習題3-4 Input: 最大負荷 2 1,110,118,8.5,20.5,0.000564,0,0,1,0,0 2,11,110,1.22,20.2,0,40,30,10,0.17,1.7 最小負荷 2 1,110,113,8.5,20.5,0.000564,0,0,1,0,0 2,11,110,1.22,20.2,0,20,15,10,0.17,1.7 南 京 理 工 大 學 《電力系統穩態分析》 課程報告 姓名 XX 學 號： 5*** 自動化學院 電氣工程 基于牛頓-拉夫遜法的潮流計算例題編程報學院(系): 專 業: 題 目: 任課教師 碩士導師 告 楊偉 XX 2015年6月10號 基于牛頓-拉夫遜法的潮流計算例題編程報告 摘要：電力系統潮流計算的目的在于：確定電力系統的運行方式、檢查系統中各元件是否過壓或者過載、為電力系統繼電保護的整定提供依據、為電力系統的穩定計算提供初值、為電力系統規劃和經濟運行提供分析的基礎。潮流計算的計算機算法包含高斯—賽德爾迭代法、牛頓-拉夫遜法和P—Q分解法等，其中牛拉法計算原理較簡單、計算過程也不復雜，而且由于人們引入泰勒級數和非線性代數方程等在算法里從而進一步提高了算法的收斂性和計算速度。同時基于MATLAB的計算機算法以雙精度類型進行數據的存儲和運算, 數據精確度高，能進行潮流計算中的各種矩陣運算，使得傳統潮流計算方法更加優化。 一 研究內容 通過一道例題來認真分析牛頓-拉夫遜法的原理和方法（采用極坐標形式的牛拉法），同時掌握潮流計算計算機算法的相關知識，能看懂并初步使用MATLAB軟件進行編程，培養自己電力系統潮流計算機算法編程能力。 例題如下：用牛頓-拉夫遜法計算下圖所示系統的潮流分布，其中系統中5為平衡節點，節點5電壓保持U=1.05為定值，其他四個節點分別為PQ節點，給定的注入功率如圖所示。計算精度要求各節點電壓修正量不大于10-6。 二 牛頓-拉夫遜法潮流計算 1 基本原理 牛頓法是取近似解x(k)之后，在這個基礎上，找到比x(k)更接近的方程的根，一步步地迭代，找到盡可能接近方程根的近似根。牛頓迭代法其最大優點是在方程f(x)=0的單根附近時誤差將呈平方減少，而且該法還可以用來求方程的重根、復根。電力系統潮流計算，一般來說，各個母線所供負荷的功率是已知的，各個節點的電壓是未知的（平衡節點外）可以根據網絡結構形成節點導納矩陣，然后由節點導納矩陣列寫功率方程，由于功率方程里功率是已知的，電壓的幅值和相角是未知的，這樣潮流計算的問題就轉化為求解非線性方程組的問題了。為了便于用迭代法解方程組，需要將上述功率方程改寫成功率平衡方程，并對功率平衡方程求偏導，得出對應的雅可比矩陣，給未知節點賦電壓初值，將初值帶入功率平衡方程，得到功率不平衡量，這樣由功率不平衡量、雅可比矩陣、節點電壓不平衡量（未知的）構成了誤差方程，解誤差方程，得到節點電壓不平衡量，節點電壓加上節點電壓不平衡量構成節點電壓新的初值，將新的初值帶入原來的功率平衡方程，并重新形成雅可比矩陣，然后計算新的電壓不平衡量，這樣不斷迭代，不斷修正，一般迭代三到五次就能收斂。2 基本步驟和設計流程圖 形成了雅克比矩陣并建立了修正方程式，運用牛頓-拉夫遜法計算潮流的核心問題已經解決，已有可能列出基本計算步驟并編制流程圖。由課本總結基本步驟如下： 1)形成節點導納矩陣Y； 2)設各節點電壓的初值，如果是直角坐標的話設電壓的實部e和虛部f；如果是極坐標的話則設電壓的幅值U和相角a； 3)將各個節點電壓的初值代入公式求修正方程中的不平衡量以及修正方程的系數矩陣的雅克比矩陣； 4)解修正方程式，求各節點電壓的變化量，即修正量； 5)計算各個節點電壓的新值，即修正后的值； 6)利用新值從第（3）步開始進入下一次迭代，直至達到精度退出循環； 7)計算平衡節點的功率和線路功率，輸出最后計算結果； ① 公式推導 ② 流程圖 三 matlab編程代碼 clear; % 如圖所示1，2，3，4為PQ節點,5為平衡節點 y=0; % 輸入原始數據，求節點導納矩陣 y(1,2)=1/(0.07+0.21j); y(4,5)=0;y(1,3)=1/(0.06+0.18j); y(1,4)=1/(0.05+0.10j); y(1,5)=1/(0.04+0.12j); y(2,3)=1/(0.05+0.10j); y(2,5)=1/(0.08+0.24j); y(3,4)=1/(0.06+0.18j); for i=1:5 for j=i:5 y(j,i)=y(i,j); end end Y=0; % 求節點導納矩陣中互導納 for i=1:5 for j=1:5 if i~=j Y(i,j)=-y(i,j); end end end % 求節點導納矩陣中自導納 for i=1:5 Y(i,i)=sum(y(i,:)); end Y % Y為導納矩陣 G=real(Y); B=imag(Y);% 輸入原始節點的給定注入功率 S(1)=0.3+0.3j; S(2)=-0.5-0.15j; S(3)=-0.6-0.25j; S(4)=-0.7-0.2j; S(5)=0; P=real(S); Q=imag(S); % 賦初值，U為節點電壓的幅值，a為節點電壓的相位角 U=ones(1,5); U(5)=1.05; a=zeros(1,5); x1=ones(8,1); x2=ones(8,1); k=0; while max(x2)>1e-6 for i=1:4 for j=1:4 H(i,j)=0; N(i,j)=0; M(i,j)=0; L(i,j)=0; oP(i)=0; oQ(i)=0; end end % 求有功、無功功率不平衡量 for i=1:4 for j=1:5 oP(i)=oP(i)-U(i)*U(j)*(G(i,j)*cos(a(i)-a(j))+B(i,j)*sin(a(i)-a(j))); oQ(i)=oQ(i)-U(i)*U(j)*(G(i,j)*sin(a(i)-a(j))-B(i,j)*cos(a(i)-a(j))); end oP(i)=oP(i)+P(i); oQ(i)=oQ(i)+Q(i); end x2=[oP,oQ]'; % x2為不平衡量列向量 % 求雅克比矩陣 % 當i~=j時，求H,N,M,L for i=1:4 for j=1:4 if i~=j H(i,j)=-U(i)*U(j)*(G(i,j)*sin(a(i)-a(j))-B(i,j)*cos(a(i)-a(j))); N(i,j)=-U(i)*U(j)*(G(i,j)*cos(a(i)-a(j))+B(i,j)*sin(a(i)-a(j))); L(i,j)=H(i,j); M(i,j)=-N(i,j); end end end % 當i=j時，求H,N,M,L for i=1:4 for j=1:5 if i~=j H(i,i)=H(i,i)+U(i)*U(j)*(G(i,j)*sin(a(i)-a(j))-B(i,j)*cos(a(i)-a(j)));N(i,i)=N(i,i)-U(i)*U(j)*(G(i,j)*cos(a(i)-a(j))+B(i,j)*sin(a(i)-a(j))); M(i,i)=M(i,i)-U(i)*U(j)*(G(i,j)*cos(a(i)-a(j))+B(i,j)*sin(a(i)-a(j))); L(i,i)=L(i,i)-U(i)*U(j)*(G(i,j)*sin(a(i)-a(j))-B(i,j)*cos(a(i)-a(j))) end end N(i,i)=N(i,i)-2*(U(i))^2*G(i,i); L(i,i)=L(i,i)+2*(U(i))^2*B(i,i); end J=[H,N;M,L] % J為雅克比矩陣 x1=-((inv(J))*x2); % x1為所求△x的列向量 % 求節點電壓新值，準備下一次迭代 for i=1:4 oa(i)=x1(i); oU(i)=x1(i+4)*U(i); end for i=1:4 a(i)=a(i)+oa(i); U(i)=U(i)+oU(i); end k=k+1; end k,U,a % 求節點注入功率 i=5; for j=1:5 P(i)=U(i)*U(j)*(G(i,j)*cos(a(i)-a(j))+B(i,j)*sin(a(i)-a(j)))+P(i); Q(i)=U(i)*U(j)*(G(i,j)*sin(a(i)-a(j))-B(i,j)*cos(a(i)-a(j)))+Q(i); end S(5)=P(5)+Q(5)*sqrt(-1); S % 求節點注入電流 I=Y*U' 四 運行結果 節點導納矩陣 經過五次迭代后的雅克比矩陣 迭代次數以及節點電壓的幅值和相角（弧度數） 節點注入功率和電流 五 結果分析 在這次學習和實際操作過程里：首先，對電力系統分析中潮流計算的部分特別是潮流計算的計算機算法中的牛頓-拉夫遜法進行深入的研讀，弄明白了其原理、計算過程、公式推導以及設計流程。牛頓-拉夫遜法是求解非線性方程的迭代過程，其計算公式為?F?J?X，式中J為所求函數的雅可比矩陣；?X為需要求的修正值；?F為不平衡的列向量。利用x(*)=x(k+1)+?X(k+1)進行多次迭代，通過迭代判據得到所需要的精度值即準確值x(*)。六 結論 通過這個任務，自己在matlab編程，潮流計算，word文檔的編輯功能等方面均有提高，但也暴漏出一些問題：理論知識儲備不足，對matlab的性能和特點還不能有一個全面的把握，對word軟件也不是很熟練，相信通過以后的學習能彌補這些不足，達到一個新的層次。 電力系統潮流計算程序設計 姓名：韋應順 學號：2011021052 電力工程學院 牛頓—拉夫遜潮流計算方法具有能夠將非線性方程線性化的特點，而使用MATLAB語言是由于MATLAB語言的數學邏輯強，易編譯。 【】【】1.MATLAB程序12 Function tisco %這是一個電力系統潮流計算的程序 n=input（‘n請輸入節點數：n=’）； m=input(‘請輸入支路數：m=’);ph=input(‘n請輸入平衡母線的節點號：ph=’)； B1=input（‘n請輸入支路信號：B1=’）;%它以矩陣形式存貯支路的情況，每行存貯一條支路 %第一列存貯支路的一個端點 %第二列存貯支路的另一個端點 %第三列存貯支路阻抗 %第四列存貯支路的對地導納 %第五列存貯變壓器的變比，注意支路為1 %第六列存貯支路的序號 B2=input(‘n請輸入節點信息：B2=’)； %第一列為電源側的功率 %第二列為負荷側的功率 %第三列為該點的電壓值 %第四列為該點的類型：1為PQ，2為PV節點，3為平衡節點 A=input（‘n請輸入節點號及對地阻抗：A=’）； ip=input（‘n請輸入修正值：ip=’）； %ip為修正值);Y=zeros（n）； Y(p,q)=Y(p,q)-1./(B1(i3)*B1(i5);e=zeros（1，n）； Y(p,q)=Y(p,q);f=zeros（1，n）； no=2*ph=1； Y(q,q)=Y(q,q)+1./B1(i3)+B1(i4)/2； End for i=1：n G=real(Y);if A（i2）=0 B=imag(Y);p=A(i1)； Y(p p)=1./A(i2);for i=1：n End e(i)=real(B2(i3));End f(i)=imag(B2(i3));For i=1:m S(i)=B2(i1)-B2(i2);p=B1(i1);V(i)=B2(i3);p=B1(i2);end Y(p,p)=Y(p,p)+1./(B1(i3)*B1(i5)^2+B1(i4)./2P=real(S);Q=imag(S);[C,D,DF]=xxf(G,B,e,f,P,Q,n,B2,ph,V,no);J=jacci(Y,G,B,P,Q,e,f,V,C,D,B2,n,ph,no);[De,Di]=hxf(J,D,F,ph,n,no);t=0;while max(abs(De))>ip&max(abs(Dfi)>ip t=t+1; e=e+De； f=f+Df; [C,D,DF]=xxf(G,B,e,f,P,Q,n,B2,ph,V,no)； J=jacci(Y,G,B,P,Q,e,f,V,C,D,B2,n,ph,no)； [De,Df]=hxf(J,Df,ph,n,no);end v=e+f*j;for i=1:n hh(i)=conj(Y(ph,i)*v(i));end S(ph)=sum(hh)*v(ph);B2(ph,1)=S(ph);V=abs(v); jd=angle(v)*180/p;resulte1=[A(:,1),real(v),imag(v),V,jd,real(S’)，imag(S’),real(B2(:1)),imag(B2(:1)),real(B2(:2)),imag(B2(:,2))];for i=1:m a(i)=conj((v(B1(i1))/B1(i5)-v(B1(i2))/B1(i3)); b(i)=v(B1(i1))*a(i)-j*B1(i4)*v(B1(i))^2/2; c(i)=-v(B1(i2))*a(i)-j*B1(i4)*v(B1(i2))^2/2;end result2=[B1(:,6),B1(:,1),B1(:,2),real(b’),imag(b’),real(c’),imag(c’), real(b’+c’),imag(b’+c’)];printcut(result1,S,b,c,result2);type resultm function [C,D,Df]=xxf(G,B,e,f,P,Q,n,B2,ph,V,no)%該子程序是用來求取Df for i=1:n If i=ph C(i)=0; D(i)=0; For j=i:n C(i)=C(i)+G(i,j)*e(j)-B(i,j)*f(j);D(i)=D(i)+G(i,j)*f(j)+B(i,j)*e(j);end P1=C(i)*e(i)+D(i)*f(i);Q1=C(i)*f(i)-D(i)*e(i);V1=e(i)^2+f(i)^2;If B2(i4)=2 p=2*i-1; Df(p)=P(i)-P1;p=p+1;else p=2*i-1; Df(p)=P(i)-P1;p=p+1; Df(p)=Q(i)-Q1;end end end Df=Df’;If ph=n Df(no?=[];end function [De,Df]=hxf(J,Df,ph,n,no)%該子函數是為求取De Df DX=JDf;DX1=DX; x1=length(DX1);if ph=n DX(no)=0;DX(no+1)=0; For i=(no+2):(x1+2)DX(i)=DX1(i-2);End Else DX=[DX1,0,0];End k=0; [x,y]=size(DX);For i=1:2:x K=k+1; Df(k)=DX(i);De(k)=DX(i+1);End End case 2 Function for j=1:n J=jacci(Y,G,B,PQ,e,f,V,C,D,B2,n,ph,no)X1=G(i,j)*f(i)-B(i,j)*e(i); X2=G(i,j)*e(i)+B(i,j)*f(i);%該子程序是用來求取jacci矩陣 for i=1:n X3=0;switch B2(i4)X4=0;case 3 P=2*i-1;continue q=2*j-1;case 1 J(p,q)=X1;for j=1:n m=p+1;if J=&J=ph J(m,q)=X3;X1=G(i)*f(i)-B(i,j)*e(i);q=q+1;X2=G(i,j)*e(i)+B(i,j)*f(i);J(p,q)=X2;X3=-X2;J(m,q)=X4;X4=X1;X1=D(i)+G(i,j)*f(i)-B(i,j)*e(i);p=2*i-1;X2=C(i)+G(i,j)*e(i)+B(i,j)*f(i);q=2*j-1;X3=0;J(p,q)=X1;X4=0;m=p+1;P=2*i-1;J(p,q)=X2;q=2*j-1;J(m,q)=X4;J(p,q)=X1;Else if j=&j=jph m=p+1;X1=D(i)+G(i,j)*f(i)-B(i,j)*e(i);J(m,q)=X3;X2=C(i)+G(i,j)*e(i)+B(i,j)*f(i);q=q+1;X3= C(i)+G(i,j)*e(i)-B(i,j)*f(i);J(p,q)=X2;X4= C(i)+G(i,j)*f(i)-B(i,j)*e(i);J(m,q)=X4;P=2*i-1;end q=2*j-1;end J(p,q)=X1;end m=p+1;end J(m,q)=X3;if ph=n q=q+1;J(no:)=[];J(p,q)=X2;J(no:)=[];J(m,q)=X4;J(:,no)=[];End J(:,no)=[];End 2實例驗證 【例題】設有一系統網絡結線見圖1，各支路阻抗和各節點功率均已以標幺值標示于圖1中，其中節點2連接的是發電廠，設節點1電壓保持U1＝1.06定值，試計算其中的潮流分布，請輸入節點數：n=5 請輸入支路數：m=7 請輸入平衡母線的節點號：ph=l 請輸入支路信息： BI=[ l 2 0．02+0．06i O l 1；1 3 0．08+0．24i 0 1 2；2 3 0．06+0．18i 0 l 3： 2 4 0．06+0．18i O l 4： 2 5 0．04+0．12i 0 l 5： 3 4 0．01+0．03i 0 l 6： 4 5 0．08+0．24i O 1 7] 請輸入節點信息： B2=[ 0 0 1．06 3；0．2+0．20i 0 1 1；一O．45一O．15i 0 l l；一0．4-0．05i 0 l 1；一0．6—0．1i 0 1 l] 請輸入節點號及對地阻抗： A=[l 0；2 0；3 0；4 0；5 O ] 請輸入修正值：ip=0．000 0l 參考文獻 [1]陳珩．電力系統穩定分析[M]．北京：中國電力出版社，2002：139—187． [2]鄭阿奇．MATLAB實用教程[M]．北京：電子工業出版社，2005：1-243． [3] 束洪春，孫士云，等．云電送粵交商流混聯系統全過 程動態電壓研究[J】．中國電力，2008，4l(10)：l-4． SHU Hong—ch吼，SUN Shi-yun，et a1．Research on fun prc'cess dyn鋤ic Voltage stabil時of hybrid AC／DC poWer tmnsmission System舶m Yu衄an proVince to G啪gdong province【J】．Electric Power，2008，4l(10)： l-4． [4] 朱新立，湯涌，等．大電網安全分析的全過程動態仿 真技術[J】．電網技術，2008，32(22)：23—28． SONG Xin—Ii，TANG Yof唱，et a1． Full dyn鋤ic simulation for the stabilhy a眥lysis of large power system【J】．Power System融IlrIolo影，2008，32(22)： 23．28． [5]Roytelm鋤I，Shallidehpour S M．A comprehcnsivc long teml dynaIIlic simulation for powcr system recoVery【J】． IEEE Transactions 0n Power Systems，1994，9(3)． [6] 石雩梅，汪志宏，等．發電機勵磁系統數學模型及參 數對電網動態穩定性分析結果影響的研究[J】．繼電 器，2007，35(21)：22-27． SHI Xue．mei，WANG Zlli-hon舀et a1．Iksearch on the innuence of g鋤e翰to璐baScd ∞de詛iled excitation system models柚d parameterS t0 power鏟id dyn鋤ic stabil時【J】．Relay，2007，35(2 1)：22-27． [7] 方思立，朱方．快速勵磁系統對系統穩定的影響[J】．中 國電機工程學報，1986，6(1)：20．28． FANG Si．1i，ZHU Fang．The effbct of f弧t．respon∞ excitation system on the stability of power netwofk【J】． Proceedings ofthe CSEE，1986，6(1)：20-28． [8] 劉取．電力系統穩定性及發電機勵磁控制[M】．北京： 中國電力出版社，2007． LIU Qu．Power system S詛bility鋤d generator excitation control【M】．BeUing：ChiIla Electric Powef Press，2007． [9] Dallachy J L，Anderson T．EXperience with rcplacing ro詛ting exciters wim static exciters【J】．1k InStitution of Electrical Engineers，1 996． [10] 陳利芳，陳天祿．淺談自并勵勵磁系統在大容量機組 中的應用【J】．繼電器，2007，35(1)：8l培4． CHEN Li-f抽島CHEN Tian—lIL Application of 辯l仁exci組tion mode in large capacity髫memtor unit【J】． ReIay'2007，35(1)：81-84． [11] 方思立，劉增煌，孟慶和．大型汽輪發電機自并勵勵 磁系統的應用條件【J]．中國電力，1994，27(12)：61．63． FANG Si．Ii，LIU Zeng-hu鋤g，MENG Qin爭hc．m application conditions of large turbine generator self-excitation system【J】．Electric Powef，1994，27(12)： 61．63． [12]梁小冰，黃方能．利用EMTDC進行長持續時間過程 的仿真研究【J】．電網技術，2002，26(9)：55．57． LIANG Xiao-bing，HUANG Fan爭眥ng．How to cany out simulalion of long dul‘ation processes by use of EMTDC【J】．Power System 11echnology，2002，26(9)： 55-57． [13]王卉，陳楷，彭哲，等．數字仿真技術在電力系統中 的應用及常用的幾種數字仿真工具【J】．繼電器，2004，32(21)：7l一75． wANG Hui，CHEN Kai，PENG zhe，et a1．Application of digital simulation眥hniques棚d severaJ simulation tools in power system[J】．Relay，2004，32(21)：71·75． [14]IEEE Power Engmeering Socie哆．IEEE std 421．5．2005 IEEE玎ccOmmended practice for excitation system models for power system stabiI時studies【s】． 電力系統潮流計算c語言程序，兩行，大家可以看看，仔細研究，然后在這個基礎上修改。謝謝 #include “stdafx.h” #include #include“Complex.h” #include“wanjing.h” #include“gauss.h” using namespace std; int _tmain(int argc, _TCHAR* argv[]){ int i; //i作為整個程序的循環變量 int N=Bus::ScanfBusNo();//輸入節點個數 int L=Line::ScanflineNo();//輸入支路個數 if((L&&N)==0){return 0;} //如果找不到兩個文件中的任意一個，退出 Line *line=new Line[L];//動態分配支路結構體 Line::ScanfLineData(line);//輸入支路參數 Line::PrintfLineData(line,L);//輸出支路參數 Bus *bus=new Bus[N];//動態分配結點結構體 for(int i=0;i bus[i].Sdelta.real=0; bus[i].Sdelta.image=0;} Bus::ScanfBusData(bus);//輸入節點參數 Bus::PrintfBusData(bus,N);//輸出結點參數 Complex **X;X=new Complex *[N];for(i=0;i Bus::JisuanNodeDnz(X,line,bus,L,N);//計算節點導納矩陣 Bus::PrintfNodeDnz(X,N);//輸出節點導納矩陣 int NN=(N-1)*2;double **JacAug;JacAug=new double *[NN];for(i=0;i double *x;x=new double[NN];int count=1; LOOP: Bus::JisuanNodeI(X,bus,N);//計算節點注入電流 Bus::JisuanNodeScal(X,bus,N);//計算節點功率 Bus::JisuanNodeScal(X,bus,N);//計算節點功率 Bus::JisuanNodeSdelta(bus,N);//計算節點功率差值 Bus::PrintfNodeScal(X,bus,N);//輸出節點功率差值 int icon=wehcon1(bus,N);//whether converbence看迭代是否結束 if(icon==1){ cout<<“icon=”< Bus::JisuanJacAug(JacAug,X,bus,N);//計算雅可比增廣矩陣 // Bus::PrintfJacAug(JacAug,N); gauss::gauss_slove(JacAug,x,NN);//解方程組求出電壓差值 Bus::ReviseNodeV(bus,x,N);//修正節點電壓 // Bus::PrintfNodeV(bus,N); count++; goto LOOP;} else { for(i=0;i { int statemp,endtemp; Complex aa,bb,cc,dd,B; B.real=0; B.image=-line[i].B; statemp=line[i].start; endtemp=line[i].end; aa=Complex::productComplex(Complex::getconj(bus[statemp-1].V), B); bb=Complex::subComplex (Complex::getconj(bus[statemp-1].V), Complex::getconj(bus[endtemp-1].V)); cc=Complex::productComplex(bb , Complex::getconj(line[i].Y)); dd=Complex::CaddC(aa,cc); line[i].stoe=Complex::productComplex(bus[statemp-1].V,dd); aa=Complex::productComplex(Complex::getconj(bus[endtemp-1].V), B); bb=Complex::subComplex (Complex::getconj(bus[endtemp-1].V), Complex::getconj(bus[statemp-1].V)); cc=Complex::productComplex(bb , Complex::getconj(line[i].Y)); dd=Complex::CaddC(aa,cc); line[i].etos=Complex::productComplex(bus[endtemp-1].V,dd); } cout<<“icon=”< Bus::JisuanNodeScal(X,bus,N);//計算節點功率 for(i=0;i { bus[i].Scal.real = bus[i].Scal.real + bus[i].Load.real;//發電機功率=注入功率+負荷功率 bus[i].Scal.image= bus[i].Scal.image+ bus[i].Load.image; bus[i].V=Complex::Rec2Polar(bus[i].V); } cout<<“====節點電壓===============發電機發出功率======”< for(i=0;i { cout<<“節點”<<(i+1)<<'t'; Complex::PrintfComplex(bus[i].V); coutt(bus[i].Scal.real); coutt(bus[i].Scal.image); cout< } cout<<“======線路傳輸功率==========”< for(i=0;i { int statemp,endtemp; statemp=line[i].start; endtemp=line[i].end; cout< Complex::PrintfComplex(Complex::ComDivRea(line[i].stoe,0.01)); Complex::PrintfComplex(Complex::ComDivRea(line[i].etos,0.01)); cout< } } return 0;} #include “stdafx.h” #include #include“Complex.h” #include“wanjing.h” #include“gauss.h” using namespace std; int _tmain(int argc, _TCHAR* argv[]){ int i; //i作為整個程序的循環變量 int N=Bus::ScanfBusNo();//輸入節點個數 int L=Line::ScanflineNo();//輸入支路個數 if((L&&N)==0){return 0;} //如果找不到兩個文件中的任意一個，退出 Line *line=new Line[L];//動態分配支路結構體 Line::ScanfLineData(line);//輸入支路參數 Line::PrintfLineData(line,L);//輸出支路參數 Bus *bus=new Bus[N];//動態分配結點結構體 for(int i=0;i bus[i].Sdelta.real=0; bus[i].Sdelta.image=0;} Bus::ScanfBusData(bus);//輸入節點參數 Bus::PrintfBusData(bus,N);//輸出結點參數 Complex **X;X=new Complex *[N];for(i=0;i Bus::JisuanNodeDnz(X,line,bus,L,N);//計算節點導納矩陣 Bus::PrintfNodeDnz(X,N);//輸出節點導納矩陣 int NN=(N-1)*2;double **JacAug;JacAug=new double *[NN];for(i=0;i double *x;x=new double[NN];int count=1; LOOP: Bus::JisuanNodeI(X,bus,N);//計算節點注入電流 Bus::JisuanNodeScal(X,bus,N);//計算節點功率 Bus::JisuanNodeScal(X,bus,N);//計算節點功率 Bus::JisuanNodeSdelta(bus,N);//計算節點功率差值 Bus::PrintfNodeScal(X,bus,N);//輸出節點功率差值 int icon=wehcon1(bus,N);//whether converbence看迭代是否結束 if(icon==1){ cout<<“icon=”< Bus::JisuanJacAug(JacAug,X,bus,N);//計算雅可比增廣矩陣 // Bus::PrintfJacAug(JacAug,N); gauss::gauss_slove(JacAug,x,NN);//解方程組求出電壓差值 Bus::ReviseNodeV(bus,x,N);//修正節點電壓 // Bus::PrintfNodeV(bus,N); count++; goto LOOP;} else { for(i=0;i { int statemp,endtemp; Complex aa,bb,cc,dd,B; B.real=0; B.image=-line[i].B; statemp=line[i].start; endtemp=line[i].end; aa=Complex::productComplex(Complex::getconj(bus[statemp-1].V), B); bb=Complex::subComplex (Complex::getconj(bus[statemp-1].V), Complex::getconj(bus[endtemp-1].V)); cc=Complex::productComplex(bb , Complex::getconj(line[i].Y)); dd=Complex::CaddC(aa,cc); line[i].stoe=Complex::productComplex(bus[statemp-1].V,dd); aa=Complex::productComplex(Complex::getconj(bus[endtemp-1].V), B); bb=Complex::subComplex (Complex::getconj(bus[endtemp-1].V), Complex::getconj(bus[statemp-1].V)); cc=Complex::productComplex(bb , Complex::getconj(line[i].Y)); dd=Complex::CaddC(aa,cc); line[i].etos=Complex::productComplex(bus[endtemp-1].V,dd); } cout<<“icon=”< Bus::JisuanNodeScal(X,bus,N);//計算節點功率 for(i=0;i { bus[i].Scal.real = bus[i].Scal.real + bus[i].Load.real;//發電機功率=注入功率+負荷功率 bus[i].Scal.image= bus[i].Scal.image+ bus[i].Load.image; bus[i].V=Complex::Rec2Polar(bus[i].V); } cout<<“====節點電壓===============發電機發出功率======”< for(i=0;i { cout<<“節點”<<(i+1)<<'t'; Complex::PrintfComplex(bus[i].V); coutt(bus[i].Scal.real); coutt(bus[i].Scal.image); cout< } cout<<“======線路傳輸功率==========”< for(i=0;i { int statemp,endtemp; statemp=line[i].start; endtemp=line[i].end; cout< Complex::PrintfComplex(Complex::ComDivRea(line[i].stoe,0.01)); Complex::PrintfComplex(Complex::ComDivRea(line[i].etos,0.01)); cout< } } return 0;} #include class Complex//定義復數類 { public: double real;double image;int RecPolar;//0表示直角坐標，1表示極坐標 static Complex CaddC(Complex c1,Complex c2);//求兩個復數和 static Complex subComplex(Complex c1,Complex c2);//求兩個復數差 static Complex productComplex(Complex c1,Complex c2);//求兩個復數積 static Complex divideComplex(Complex c1,Complex c2);//求兩個復數商 static Complex ComDivRea(Complex c1,double r2);//除數 static Complex getconj(Complex c1);//求一個復數共軛 static Complex getinverse(Complex c1);//取倒數 static double getComplexReal(Complex c1);//求一個復數實部 static double getCompleximage(Complex c1);//求一個復數虛部 static void PrintfComplex(Complex c1);//顯示一個復數 static void PrintfmultiComplex(Complex C,int N);//顯示多個復數 static void zeroComplex(Complex c1);//將復數復零 static Complex Rec2Polar(Complex c1);//取極坐標 Complex(){ RecPolar=0;} }; Complex Complex::Rec2Polar(Complex c1)//極坐標表示 { Complex Node;Node.real=sqrt(c1.real*c1.real+c1.image*c1.image);Node.image=atan2(c1.image,c1.real)*180/3.1415926;Node.RecPolar=1;return Node;} Complex Complex::CaddC(Complex c1,Complex c2)//復數加法 { Complex Node; Node.real=c1.real+c2.real; Node.image=c1.image+c2.image; return Node;} Complex Complex::subComplex(Complex c1,Complex c2)//復數減法 { Complex Node; Node.real=c1.real-c2.real; Node.image=c1.image-c2.image; return Node;} Complex Complex::productComplex(Complex c1,Complex c2)//復數乘法 { Complex Node; Node.real=c1.real*c2.real-c1.image*c2.image; Node.image=c1.image*c2.real+c2.image*c1.real; return Node;} Complex Complex::divideComplex(Complex c1,Complex c2)//復數除法 { Complex Node; Node.real=(c1.real*c2.real+c1.image*c2.image)/(pow(c2.real,2)+pow(c2.image,2));Node.image=(c1.image*c2.real-c1.real*c2.image)/(pow(c2.real,2)+pow(c2.image,2));return Node;} Complex Complex::ComDivRea(Complex c1,double r1)//復數除數 { Complex Node;Node.real=c1.real/(r1);Node.image=c1.image/(r1);return Node;} Complex Complex::getconj(Complex c1)//取共軛 { Complex Node; Node.real=c1.real;Node.image=-c1.image; return Node;} Complex Complex::getinverse(Complex c1)//取倒數 { Complex Node;Node.real=1;Node.image=0;Node=(Complex::divideComplex(Node,c1));return Node;} double Complex::getComplexReal(Complex c1)//取實部 { return c1.real;} double Complex::getCompleximage(Complex c1)//取虛部 { return c1.image;} void Complex::PrintfComplex(Complex c1)//按直角坐標輸出 { if(c1.RecPolar==0){ cout.precision(6); cout.width(8); cout.setf(ios::right); cout< ”; cout.precision(6); cout.width(8); cout.setf(ios::left); cout< ”;} else { cout< Complex::zeroComplex(Complex c1)//清零 { c1.real=0;c1.image=0;} class gauss { public: static void gauss_slove(double **a,double *x,int NN);static void gauss_output();}; void gauss::gauss_slove(double **a,double *x,int NN){ int n,i,j,k,*pivrow,**pivcol;double eps,pivot,sum,aik,al; n=NN;pivrow=new int[n];pivcol=new int *[n]; for(i=0;i pivot= fabs(a[k][k]); pivrow[k]=k;//行 pivcol[k][0]=k;pivcol[k][1]=k;//列n*2矩陣 for(i=k;i { for(j=k;j { if(pivot { pivot=fabs(a[i][j]); pivrow[k]=i;//行 pivcol[k][1]=j;//列 } } } if(pivot { cout<<“error”< getchar(); exit(0); } if(pivrow[k]!=k)//行變換 { for(j=k;j<(n+1);j++) { al=a[pivrow[k]][j]; a[pivrow[k]][j]=a[k][j]; a[k][j]=al; } } if(pivcol[k][1]!=k)//列變換 { for(i=0;i { al=a[i][pivcol[k][1]]; a[i][pivcol[k][1]]=a[i][k]; a[i][k]=al; } } if(k!=(n-1))//將矩陣化為上三角形 式 { for(i=(k+1);i { aik=a[i][k]; for(j=k;j<(n+1);j++) { a[i][j]-=aik*a[k][j]/a[k][k]; } } } } x[n-1]=a[n-1][n]/a[n-1][n-1];//解方程 for(i=(n-2);i>=0;i--){ sum=0; for(j=(i+1);j { sum +=a[i][j]*x[j];0.182709 0.016894-0.0310701 -0.0402051 0.156702 -0.0355909-0.0668055 -0.00703229-0.0886481 -0.0129814-0.0390805 -0.0135062-0.1023 -0.0460568 -0.0342827 -0.00382402-0.102896 -0.0184062 } x[i]=(a[i][n]-sum)/a[i][i];} for(k=(n-2);k>=0;k--){ al=x[pivcol[k][1]]; x[pivcol[k][1]]=x[pivcol[k][0]]; x[pivcol[k][0]]=al;} cout<<“節點電壓修正量”< cout< } ====節點功率計算值==== 0.935261 -0.159048 0.573909 0.0789973-0.00289889 -0.00796623-0.0791247 -0.0168362-0.436255 -0.0580392 0.0359139 -0.0106592-0.229118 -0.0885419-0.136179 -0.148207 0.0446243 0.0111298-0.0223764 -0.00695775-0.0237482 -0.198318 -5.24266e-015 -0.0354071 -0.0925078 -1.05629e-015 -0.0391348 0.014529 0.00158644 -0.0258771 -0.109514 icon=1進行第2次迭代 節點電壓修正量 =================-0.00164889-0.000540034-0.00261067-0.00532027-0.00235315-0.00600971-0.00189677-0.00643874-0.0023631-0.00650659-0.00170949-0.0074907-0.00164545-0.00485415-0.00493977-0.0119042-0.00331285-0.0175611-0.00207908 -0.00347744-0.0869347-9.48909e-015-0.0110778-0.0538236-7.53784e-016-0.0168097 7.049e-005-0.00146487-0.00458276 0.00251645 -0.00336375-0.00530645-0.0147816-0.000326161-0.00640487-0.00251701-0.0169829-0.00175286-0.0174333-0.0239063 -0.0119192-0.076014 -0.0160104-0.441997 -0.0750285 0.000250012 3.72542e-005-0.228052 -0.108844-0.100078 -0.105634 0.000410707 0.000378067-0.057497 -0.0195879 0.200039 0.0582563-0.00307326-0.0163809-0.00232773-0.0175806 8.74293e-005-0.0192018 0.000558996-0.0197776-0.000247851-0.0193784-0.00115346-0.0185848-0.00127275-0.0186244-0.00010108-0.0188966 0.000553585-0.0200901-3.76315e-005-0.0208303 0.00308341-0.0219386-0.00195916-0.0205356-0.00184757-0.0076401 0.00197593-0.0245534 0.00434657-0.027534 ====節點功率計算值==== 0.98623 -0.134163 0.583136 0.166278-0.111173 0.199792 -0.0621041 -0.0821379 -0.0350785 -0.0902383 -0.0320461 -0.0951562 -0.0220362 -0.175458 4.72557e-015 -0.0320661 -0.0871134 -7.03489e-017 -0.0350769 0.000273455 1.51804e-005 -0.0240417 -0.10604 icon=1進行第3次迭代 節點電壓修正量 =================-2.67079e-005-2.30128e-006-2.20543e-005-6.00686e-005-2.33043e-005-6.85601e-005-3.22294e-005-2.61107e-005-2.80198e-005-6.6167e-005-2.34528e-005 -0.0739846 0.0227868-0.0158709-0.0248173-0.0179447-0.0578368-0.00890719-0.0337091-0.00693706-0.111601 1.21429e-014-0.0159145-0.0667319 9.24355e-016-0.0228592 7.10354e-005-6.6188e-006-0.00889343-0.0184098 -5.66132e-005-4.4646e-005-1.74668e-005-4.50947e-005-0.000181763-3.81763e-006-0.000286581-6.68993e-005-1.28441e-005-5.17172e-005-0.000223284-4.54717e-005-2.47586e-005 4.32335e-007-0.000258494 1.82635e-005-0.000272051-6.95195e-006-0.000251969 1.11318e-005-0.000279418 5.74737e-005-0.000307368 6.86998e-005-0.000320274 5.38112e-005-0.00031447 3.59531e-005-0.00030494 3.37607e-005-0.000307449 5.26532e-005-0.000310721 6.92761e-005-0.000350373 5.60942e-005-0.00040977 0.000123641-0.000440259 1.36149e-005-0.000426973-1.70227e-005-9.37794e-005 0.000113675-0.000544011 0.000176034-0.000636202 ====節點功率計算值==== 0.986878 -0.133979 0.583 0.167193-0.024 -0.012-0.076 -0.016-0.442 -0.0748606 1.43501e-008 1.07366e-008-0.228 -0.109 -0.0999999 -0.104049 4.51318e-008 8.98835e-008-0.0579999 -0.0199999 0.2 0.0591018-0.112 -0.0749997 0.2 0.0242519-0.062 -0.016-0.082 -0.025-0.035 -0.018 -0.0900001 -0.058-0.032 -0.00899997-0.095 -0.0339999-0.022 -0.00699998-0.175 -0.112 -6.07156e-015 -1.19217e-014-0.032 -0.016-0.087 -0.0669999 7.03078e-017 -9.23979e-016-0.035 -0.0229999 1.09492e-007 4.45699e-008 1.54958e-009 -2.01531e-010-0.024 -0.00899994-0.106 -0.0189996 icon=0,迭代結束。 ====節點電壓===============發電機發出功率====== 節點1 1.05 0。 98.6878-13.3979 節點2 1.045 -1.846。 29.4193 節點3 1.02384-3.83352。 0 節 點25 1.01216-9.68486。 0 0 0 節點4 1.01637-4.55698。 0 節 點26 0.994393 -10.1089。 0 0 0 節點5 1.01 -6.48617。 節 點27 1.02012-9.42025。 0 11.5139 0 節點6 1.01332-5.38073。 0 節 點28 1.00992-5.86244。 0 0 0 節點7 1.00489-6.38368。 0 節 點29 1.00022-10.6579。 0 0 節點8 19.5951 節點9 0 節點10 0 節點11 5.91018 節點12 0 節點13 2.42519 節點14 0 節點15 0 節點16 0 節點17 0 節點18 0 節點19 0 節點20 0 節點21 0 節點22 0 節點23 0 節點24 0 1.01 -5.62974。 1.03905-6.78143。 1.03595-8.69362。 -4.5962。 1.04711-7.80323。 1.05 -6.34392。 1.03242-8.7401。 1.02788-8.86784。 1.03458-8.45044。 1.03051-8.83678。 1.01845-9.5141。 1.01604-9.70326。 1.02022-9.50938。 1.0237-9.17478。 1.02432-9.17024。 1.01802-9.36719。 1.01339-9.68362。 0 20 節 點30 0.988705 -11.5464。 0 0 0 ====== 線路傳輸功率========== 2to1 -57.7373 5.41674i 58.3454 0 -15.1827i 3to1 -39.659 -7.75964i 40.3424 1.78481i 4to2 -30.87 -9.74186i 31.4153 0 3.58352i 4to3 -37.0772 -7.78596i 37.259 6.55964i 5to2 -44.3717 -9.78456i 45.2968 0 4.84242i 6to2 -38.4766 -8.22625i 39.3252 0 2.87667i 6to4 -34.946 1.92384i 35.0885 0 -3.28202i 7to5 -0.16304 -6.41767i 0.171702 0 2.2985i 7to6 -22.637 -4.48233i 22.7745 0 1.44238i 8to6 -11.8939 -5.48098i 11.913 0 3.70557i 6to9 12.3737 -12.3826i -12.3737 0 13.0033i 6to10 10.9107 -3.80907i -10.9107 0 4.53223i 11to9 5.91018i 0 -5.08963i 10to9 -32.652 -2.3712i 32.652 0 3.46974i 4to12 23.5411 -11.5375i -23.5411 0 13.2407i 13to12 2.42519i 1.05 -1.90978i 1.66484i 14to12 -7.9019 -2.06732i 7.97894 30to29 -3.6702 -0.542564i 3.70398 2.22749i 0.606393i 15to12 -18.254 -5.74885i 18.4835 28to8 -1.89152 -3.79982i 1.89395 6.20089i-4.9239i 16to12-7.53872 -2.90237i 7.59633 28to6 -14.7868 -2.82565i 14.8234 3.02352i 0.294601i 15to14-1.69544 -0.461488i 1.70189 請按任意鍵繼續...0.467323i 17to16-4.03014 1.10238i 18to15-6.08074 1.46028i 19to18-2.87549 0.478389i 20to19 6.6418-2.93222i 20to10 -8.8418 3.85077i 17to10-4.96987 4.76656i 21to10-16.1562 9.42843i 22to10-7.87782 4.21401i 22to21 1.34443-2.01837i 23to15-5.59369 2.25006i 24to22-6.48186 2.08163i 24to23-2.38596 0.579814i 25to24-0.167617 0.281364i 26to25 -3.5 2.3674i 27to25 3.39433-2.08638i 28to27 16.1446 3.13006i 29to27-6.10398 1.67047i 30to27-6.92979-1.07089i-1.37839i-0.467767i 2.96679i-3.66679i-4.72911i-9.18162i-4.10132i 2.01969i-2.17981i-2.00141i-0.56401i -0.28102i-2.29999i 2.11848i-2.10093i-1.50639i -1.3574i 4.03872 6.12096 2.88074 -6.62452 8.9242 4.98423 16.2709 7.93248 -1.34378 5.62846 6.53339 2.39369 0.167814 3.54513 -3.37751 -16.1446 6.19083 7.09313 高等電力系統分析 IEEE30節點潮流程序 班級：電研114班 姓名：王大偉 學號：2201100151 電力系統潮流計算發展史 對潮流計算的要求可以歸納為下面幾點： （1）算法的可靠性或收斂性（2）計算速度和內存占用量（3）計算的方便性和靈活性 電力系統潮流計算屬于穩態分析范疇，不涉及系統元件的動態特性和過渡過程。因此其數學模型不包含微分方程，是一組高階非線性方程。非線性代數方程組的解法離不開迭代，因此，潮流計算方法首先要求它是能可靠的收斂，并給出正確答案。隨著電力系統規模的不斷擴大，潮流問題的方程式階數越來越高，目前已達到幾千階甚至上萬階，對這樣規模的方程式并不是采用任何數學方法都能保證給出正確答案的。這種情況促使電力系統的研究人員不斷尋求新的更可靠的計算方法。 在用數字計算機求解電力系統潮流問題的開始階段，人們普遍采用以節點導納矩陣為基礎的高斯-賽德爾迭代法（一下簡稱導納法）。這個方法的原理比較簡單，要求的數字計算機的內存量也比較小，適應當時的電子數字計算機制作水平和電力系統理論水平，于是電力系統計算人員轉向以阻抗矩陣為主的逐次代入法（以下簡稱阻抗法）。 20世紀60年代初，數字計算機已經發展到第二代，計算機的內存和計算速度發生了很大的飛躍，從而為阻抗法的采用創造了條件。阻抗矩陣是滿矩陣，阻抗法要求計算機儲存表征系統接線和參數的阻抗矩陣。這就需要較大的內存量。而且阻抗法每迭代一次都要求順次取阻抗矩陣中的每一個元素進行計算，因此，每次迭代的計算量很大。 阻抗法改善了電力系統潮流計算問題的收斂性，解決了導納法無法解決的一些系統的潮流計算，在當時獲得了廣泛的應用，曾為我國電力系統設計、運行和研究作出了很大的貢獻。但是，阻抗法的主要缺點就是占用計算機的內存很大，每次迭代的計算量很大。當系統不斷擴大時，這些缺點就更加突出。為了克服阻抗法在內存和速度方面的缺點，后來發展了以阻抗矩陣為基礎的分塊阻抗法。這個方法把一個大系統分割為幾個小的地區系統，在計算機內只需存儲各個地區系統的阻抗矩陣及它們之間的聯絡線的阻抗，這樣不僅大幅度的節省了內存容量，同時也提高了節省速度。 克服阻抗法缺點的另一途徑是采用牛頓-拉夫遜法（以下簡稱牛頓法）。牛頓法是數學中求解非線性方程式的典型方法，有較好的收斂性。解決電力系統潮流計算問題是以導納矩陣為基礎的，因此，只要在迭代過程中盡可能保持方程式系數矩陣的稀疏性，就可以大大提高牛頓潮流程序的計算效率。自從20世紀60年代中期采用了最佳順序消去法以后，牛頓法在收斂性、內存要求、計算速度方面都超過了阻抗法，成為直到目前仍被廣泛采用的方法。 在牛頓法的基礎上，根據電力系統的特點，抓住主要矛盾，對純數學的牛頓法進行了改造，得到了P-Q分解法。P-Q分解法在計算速度方面有顯著的提高，迅速得到了推廣。 牛頓法的特點是將非線性方程線性化。20世紀70年代后期，有人提出采用更精確的模型，即將泰勒級數的高階項也包括進來，希望以此提高算法的性能，這便產生了保留非線性的潮流算法。另外，為了解決病態潮流計算，出現了將潮流計算表示為一個無約束非線性規劃問題的模型，即非線性規劃潮流算法。 近20多年來，潮流算法的研究仍然非常活躍，但是大多數研究都是圍繞改進牛頓法和P-Q分解法進行的。此外，隨著人工智能理論的發展，遺傳算法、人工神經網絡、模糊算法也逐漸被引入潮流計算。但是，到目前為止這些新的模型和算法還不能取代牛頓法和P-Q分解法的地位。由于電力系統規模的不斷擴大，對計算速度的要求不斷提高，計算機的并行計算技術也將在潮流計算中得到廣泛的應用，成為重要的研究領域。第二篇：電力系統潮流計算
第三篇：電力系統潮流計算程序設計
第四篇：電力系統潮流計算程序
第五篇：電力系統潮流計算發展史
點擊咨詢 