no2.txt

本程序采用牛顿拉夫逊进行电力网络的潮流分布计算

#define L 14				//网络的支路总数
#define N 8					//节点数
#define M 4					//网络的PQ节点数
#define N0 2*(N-1)			//雅可比矩阵的行数
#define N1 N0+1
int _tmain(int argc, _TCHAR* argv[])
{
	int K=0;					//打印开关.K=1,则打印;否则,不打印
	int K1=0;					//子程序PLSC中判断输入电压的形式;K1=1,则为极座标形式.否则为直角坐标形式. 
	float D=1e-4;				//有功及无功功率误差的最大值（绝对值）.
	float *DD;					//误差量最大者
	float p=1.0;
	DD=&p;
	float G[f2(N,N,N)+1];				//Ybus的电导元素(实部).  
	float B[f2(N,N,N)+1];				//Ybus的电纳元素(虚部).   
	float G1[L]={0.2,0.16,0.35,0.3,0.21,0.2,0.21,0.32,0.25,0.1,1.1,1.0,0.75,0.8};	//第I支路的串联电导.      
	float B1[L]={1.3,2.1,2.0,3.0,1.0,1.0,3.1,3.0,2.3,1.5,8.1,7.0,6.3,6.5};	//第I支路的串联电纳
	float C1[L]={0,0,0,0,0,0,0,0,0,0,0.08,0.05,0.06,0.03};			//C1(I)第I支路的pie型对称接地电纳
	float C[f2(N,L,L)+1]={0};	//C(I,J)第I节点J支路不对称接地电纳
	float CO[N]={0};				//CO(I) :第I节点的接地电纳
	int S1[L]={1,1,1,1,2,2,3,3,3,6,5,5,6,7};					//第I支路的起始节点号
	int E1[L]={2,4,7,8,3,8,4,5,6,7,6,8,7,8};					//第I支路的终止节点号
	float P[N]={0.25,-0.22,0.25,0,-0.233,0.15,0.20,0};				//第I节点的注入有功功率
	float Q[N]={0.2,-0.13,0,-1.0,0,0,0,0};				//第I节点的注入无功功率
	float P0[N]={0};				//第I节点有功功率误差
	
	float V0[N]={0};				//第I节点(PV节点)的电压误差(平方误差)
	float V[N]={1,1,1,1,1,1,1,1};				//第I节点的电压幅值
	
	float E[N]={1,1,1,1,1,1,1,1};				//第I节点的电压的实部
	float F[N]={0,0,0,0,0,0,0,0};				//第I节点的电压的虚部
	float JM[f2(N0,N0,N0)+1];			//Jacoby矩阵的第I行J列元素
	float A[f2(2*N-2,2*N-1,2*N-1)+1]={0};		//A(I,J)修正方程的增广矩阵,三角化矩阵的第I行J列元素,
										//运算结束后A矩阵的最后一列存放修正的解.  
	float P1[L];				//第I支路由S1(I)节点注入的有功功率
	float Q1[L];				//第I支路由S1(I)节点注入的无功功率
	float P2[L];				//第I支路由E1(I)节点注入的有功功率
	float Q2[L];				//第I支路由E1(I)节点注入的无功功率
	float P3[L];				//第I支路的有功功率损耗
	float Q3[L];				//第I支路的无功功率损耗
	float ANGLE[N];			//第I节点电压的角度
	float Q0[N]={0};			//第I节点有功功率误差
	int ii,jj;
	int cnt=1;				//迭代次数
//	file4=fopen("cs_flow_1.txt","w");
	ybus(N,L,M,G,B,G1,B1,C1,C,CO,1,S1,E1);//计算电导和导纳参数
	dpqc(P,Q,P0,Q0,V,V0,M,N,E,F,1,G,B,DD);
	//while (arrayMax(N,Q0)>D || arrayMax(N,P0)>D){
	while (*(DD)>D){
		cnt++;
		jmcc(M,N,N0,E,F,G,B,JM,1);//算jacob矩阵
		//算增广的修正方程矩阵，最后一列放误差向量，与书上的方程顺序稍有不同
		for (ii=1;ii<=N0;ii++){
			for(jj=1;jj<=N0;jj++){
				A[f2(ii,jj,N1)]=JM[f2(ii,jj,N0)];
			}
		}
		for(ii=1;ii<=N0;ii++){
			if(ii<=2*M){
				if(ii%2==1)		A[f2(ii,2*N-1,2*N-1)]=Q0[f1((ii+1)/2)];
				else	A[f2(ii,2*N-1,2*N-1)]=P0[f1(ii/2)];
			}
			else{
				if(ii%2==1)		A[f2(ii,2*N-1,2*N-1)]=V0[f1((ii+1)/2)];
				else	A[f2(ii,2*N-1,2*N-1)]=P0[f1(ii/2)];
			}
		}
		sevc(A,N0,1,N1);//解线性方程组
		for(ii=1;ii<=N0;ii++){
			if(ii%2==1){
				E[f1((ii+1)/2)]-=A[f2(ii,2*N-1,2*N-1)];
			}
			else F[f1(ii/2)]-=A[f2(ii,2*N-1,2*N-1)];
		}//修正解
		dpqc(P,Q,P0,Q0,V,V0,M,N,E,F,1,G,B,DD);//重新计算误差
	}

	/*float Atest[f2(4,5,5)+1]={1,2,0,0,1, 0,1,0,0,2, 0,0,1,0,1, 2,0,0,1,1};
	sevc ( Atest, 4, 1, 5);*/
	plsc(  N,  L,  M,  G ,  B ,  E ,  F ,\
             E1 ,  S1 ,  G1 ,  B1 ,  C1 ,  C ,\
             CO ,  P1 ,  Q1 ,  P2 ,  Q2 ,  P3 ,\
             Q3 ,  P ,  Q ,  V ,  ANGLE , K1);
	//file4=fclose("cs_flow_1");
	return 0;
}