cnetflow.cpp

电力系统潮流程序adsfasdfasfd

#include "CNetFlow.h"
#include "db_def/db_define.h"
#include "paramanage.h"
#include <assert.h>
#include <algorithm>
#include <math.h>
#include <fstream>
using namespace NM_PSNTP;
#include <fstream>
const float V0 = 1.0;
const float THETA0 = 0.0;
CNetFlow::CNetFlow()
{
	Branch_T branch;
//支路的输入排序按i的大小从小到大排序。
	
	branch.i = 0;
	branch.j = 1;
	branch.r = 0.04;
	branch.x = 0.25;
	branch.yk = 0.5;
	m_allflow_dev_vec.branch_vec.push_back(branch);	
	
	branch.i = 0;
	branch.j = 2;
	branch.r = 0.1;
	branch.x = 0.35;
	branch.yk = 0;
	m_allflow_dev_vec.branch_vec.push_back(branch);

	branch.i = 1;
	branch.j = 2;
	branch.r = 0.08;
	branch.x = 0.30;
	branch.yk = 0.5;//这里填的b0是线路整个的，没有除以2的
	m_allflow_dev_vec.branch_vec.push_back(branch);	

	branch.i = -1;
	branch.j = 3;
	branch.r = 0.0;
	branch.x = 0.015;
	branch.yk = 1.05;
	m_allflow_dev_vec.branch_vec.push_back(branch);	
	
	branch.i = -2;
	branch.j = 4;
	branch.r = 0;
	branch.x = 0.03;
	branch.yk = 1.05;
	m_allflow_dev_vec.branch_vec.push_back(branch);	
	

	SetNodeNum(5);
	m_DI.resize(m_node_num);//m_DI存储的是ax=b的b值,等于节点注入功率-节点计算功率
	
	PvNode_T pvnode;
	pvnode.i = 3;
	pvnode.v = 1.05;
	m_allflow_dev_vec.pvnode_vec.push_back(pvnode);

	pvnode.i = 4;
	pvnode.v = 1.05;
	m_allflow_dev_vec.pvnode_vec.push_back(pvnode);
	
	Generator_T gen_info;
	Load_T load_info;
	load_info.i = 0;
	load_info.p = -1.6;
	load_info.q = -0.8;
	load_info.v = 1.0;
    m_allflow_dev_vec.load_vec.push_back(load_info);
	
	load_info.i = 1;
	load_info.p = -2;
	load_info.q = -1.3;
	load_info.v = 1.0;
	m_allflow_dev_vec.load_vec.push_back(load_info);
	
	load_info.i = 2;
	load_info.p = -3.7;
	load_info.q = -1.3;
	load_info.v = 1.0;
	m_allflow_dev_vec.load_vec.push_back(load_info);
	
	gen_info.i = 3;
	gen_info.p = 5;
	gen_info.q = 0;
	gen_info.v = 1.05;
	m_allflow_dev_vec.gen_vec.push_back(gen_info);
	
	
	gen_info.i = 4;
	gen_info.p = 0;
	gen_info.q = 0;
	gen_info.v = 1.05;
	m_allflow_dev_vec.gen_vec.push_back(gen_info);
	

	FormAdmittanceMatrix();
	InitNodeVolt();
    m_isConv = false;
	m_itertimes = 0;
	while(!m_isConv)
	{
        if(m_itertimes >10)
			break;
		FormFactorTable(1);
		CalNodePower(1);
		BackSubstitution();

		FormFactorTable1(2);
		CalNodePower(2);
		BackSubstitution();

		AmendVolt(1);
		AmendVolt(2);
		CalError(1);
		CalError(2);
		m_itertimes++;
	}
	CalBranchFlow();
}       

CNetFlow::~CNetFlow()
{

}
void CNetFlow::SetNodeNum(int node_num)
{
    m_node_num = node_num;
} 
	
void CNetFlow::Reset()
{
    m_yii_vec.clear();
    m_yii1_vec.clear();
    m_yij1_vec.clear();
    m_yij_vec.clear();
    m_nysum_vec.clear();
    m_nyseq_vec.clear();

    m_node_num = 0;//系统中"母线"个数，也就是N
    m_branch_num = 0;//系统中"支路"个数，也就是线路+变压器
    m_gn_num = 0;//发电机个数
    m_ld_num = 0;//负荷个数
  
    m_node_p.clear();;
    m_node_q.clear();;
    m_ndvolt_vec.clear();;
  
    m_DI.clear();;
    m_nusum_vec.clear();;//存储因子表中的上三角矩阵的每行非零元素个数
    m_d_vec.clear();//存储因子表对角元素
    m_u_vec.clear();//存储因子表上三角非零元素。
}

void CNetFlow::SetGraph(AllFlowDevice_T &alldev)
{
 	m_allflow_dev_vec = alldev;
}
//形成节点导纳矩阵
void CNetFlow::FormAdmittanceMatrix()
{
	int         i,j,n;
	float       r,x,yk;
	float       zmag2,gij,bij;
	float       b_ij;
	m_yii_vec.resize(m_node_num);
	m_yii1_vec.resize(m_node_num);
	m_nysum_vec.resize(m_node_num);
	m_nyseq_vec.resize(m_node_num);
	//先初始化为0
	for(i = 0;i<m_node_num;++i)
	{
		m_yii_vec[i].g = 0;
		m_yii_vec[i].b = 0;
		m_yii1_vec[i].g = 0;
		m_yii1_vec[i].b = 0;
		m_nysum_vec[i] = 0;
	}
	
	m_branch_num = m_allflow_dev_vec.branch_vec.size();
	m_yij_vec.resize(m_branch_num);
	m_yij1_vec.resize(m_branch_num);
	//不考虑接地支路，形成节点导纳矩阵
	for(n = 0;n < m_branch_num;++n)
	{
		i = abs(m_allflow_dev_vec.branch_vec[n].i);
		j = abs(m_allflow_dev_vec.branch_vec[n].j);
		r = m_allflow_dev_vec.branch_vec[n].r;
		x = m_allflow_dev_vec.branch_vec[n].x;
		yk = m_allflow_dev_vec.branch_vec[n].yk;
		
		zmag2 = r*r+x*x;
		gij = r/zmag2;
		bij = -x/zmag2;//书本上这个地方缺个负号　
		
		b_ij = -1.0/x;
		
		//如果是变压器支路
		if(m_allflow_dev_vec.branch_vec[n].i< 0 || m_allflow_dev_vec.branch_vec[n].j <0)
		{
			m_yij_vec[n].g = -gij/yk;
			m_yij_vec[n].b = -bij/yk;
			m_yij1_vec[n].g = 0.0;
			m_yij1_vec[n].b = -b_ij/yk;
		}
		//如果是线路支路
		else
		{
			m_yij_vec[n].g = -gij;
			m_yij_vec[n].b = -bij;
			m_yij1_vec[n].g = 0.0;
			m_yij1_vec[n].b = -b_ij;		
		}
		m_yij_vec[n].j = j;
		m_yij1_vec[n].j = j;
	
		//如果是变压器支路
		if(m_allflow_dev_vec.branch_vec[n].i< 0 || m_allflow_dev_vec.branch_vec[n].j <0)
		{
			m_yii_vec[i].g -= m_yij_vec[n].g;
			m_yii_vec[i].b -= m_yij_vec[n].b ;
			m_yii_vec[j].g -= m_yij_vec[n].g;
			m_yii_vec[j].b -= m_yij1_vec[n].b;
			
		    m_yii1_vec[i].b += b_ij/yk;
			m_yii1_vec[j].b += b_ij/yk;
		}
		//如果是线路支路
		else
		{
			m_yii_vec[i].g += gij;
			m_yii_vec[i].b += bij;
			m_yii_vec[j].g += gij;
			m_yii_vec[j].b += bij;
			
		    m_yii1_vec[i].b += b_ij;
			m_yii1_vec[j].b += b_ij;
		}
		m_nysum_vec[i]++;
	}
	
	m_nyseq_vec[0] = 0;//书本上是从1开始，我的数组是从0开始。
	for(i = 0;i < m_node_num-1;++i)
	{
		m_nyseq_vec[i+1] = m_nyseq_vec[i]+m_nysum_vec[i];
		printf("i====%d,m_nyseq_vec[i+1]===%d\n",i+1,m_nyseq_vec[i+1]);
	}
	printf("m_yii_vec[i].g====%f,m_yii_vec[i].b===%f\n",m_yii_vec[0].g,m_yii_vec[0].b);
	
	//下面在导纳矩阵中加入接地之路
	for(n = 0;n < m_branch_num;++n)
    {
  	    i = m_allflow_dev_vec.branch_vec[n].i;
  	    j = m_allflow_dev_vec.branch_vec[n].j;
  	    yk = m_allflow_dev_vec.branch_vec[n].yk;
  	
  	    //变压器支路
  	    if(i < 0 || j < 0)
  	    {
  		if(i < 0)
  		{
  			i = abs(i);
  			gij = m_yij_vec[n].g;
  			bij = m_yij_vec[n].b;
  			
  			b_ij =  m_yij1_vec[n].b;
  			
  			m_yii_vec[i].g += (1.0 -1.0/yk)*gij;
  			m_yii_vec[i].b += (1.0 -1.0/yk)*bij;
  			
  			m_yii1_vec[i].b += (1.0 -1.0/yk)*b_ij;
  			
  			m_yii_vec[j].g += (1.0 -yk)*gij;
  			m_yii_vec[j].b += (1.0 -yk)*bij;
  			
  			m_yii1_vec[j].b += (1.0 -yk)*b_ij;
  		}
  		else
  		{
  			j = abs(j);
  			gij = m_yij_vec[n].g;
  			bij = m_yij_vec[n].b;
  			
  			b_ij =  m_yij1_vec[n].b;
  			
  			m_yii_vec[j].g += (1.0 -1.0/yk)*gij;
  			m_yii_vec[j].b += (1.0 -1.0/yk)*bij;
  			
  			m_yii1_vec[j].b += (1.0 -1.0/yk)*b_ij;
  			
  			m_yii_vec[i].g += (1.0 -yk)*gij;
  			m_yii_vec[i].b += (1.0 -yk)*bij;
  			
  			m_yii1_vec[i].b += (1.0 -yk)*b_ij;
  		}
  	}
  	//线路支路
  	else
  	{
  		bij = yk/2.0;
  		b_ij = yk/2.0;
  		m_yii_vec[i].b += bij;
   		m_yii_vec[j].b += bij;
   		
   		m_yii1_vec[i].b += b_ij;
   		m_yii1_vec[j].b += b_ij;  		
  	}
  }
  for(i = 0;i<m_yii_vec.size();++i)
  {
	  printf("i == %d,g===%f,b===%f\n",i,m_yii_vec[i].g,m_yii_vec[i].b);
  }
  for(i = 0;i<m_yij_vec.size();++i)
  {
	  printf("j == %d,g===%f,b===%f\n",m_yij_vec[i].j,m_yij_vec[i].g,m_yij_vec[i].b);
  }
}

//形成因子表，flag = 1形成B',flag = 2,形成B"
void CNetFlow::FormFactorTable(int flag)
{
	 int n_pv = 0;
	 int n_u;//因子表上三角矩阵元素计数变量
	 int count,i_above;//消去行号的计数变量
	 int i_pv = m_allflow_dev_vec.pvnode_vec[0].i;
	 int i,j;
	 vector<float> B;
	 float Btemp;
	 B.resize(m_node_num);
	 m_u_vec.resize(m_branch_num);//???
     m_d_vec.resize(m_node_num);
	 m_nusum_vec.resize(m_node_num);
	 for(i = 0;i < m_node_num-1;++i)
	 {
	 	  //如果是计算B",并且是pv节点
	 	  if(2 == flag && i == i_pv)
	 	  {
	 	  	n_pv ++;
	 	  	i_pv = m_allflow_dev_vec.pvnode_vec[count].i;
	 	  	m_nusum_vec[i] = 0;
	 	  	m_d_vec[i] = 0;
	 	  	
	 	  }
	 	  else
	 	  {
	 	  	for(count = i+1;count < m_node_num-1;++count)
	 	  	{
	 	  		 B[count]= 0.0;	 
	 	  	}
	 	  	B[i] = m_yii_vec[i].b;
			float temp1 = B[i];
	 	  	for(count  = m_nyseq_vec[i];count <m_nyseq_vec[i+1];++count)
	 	  	{
	 	  		j = m_yij_vec[count].j;
	 	  		B[j] = m_yij_vec[count].b;
				float temp1 = B[j];
	 	  	}
	 	  	if(2 == flag)
	 	  	{
	 	  		for(count  = 0;count < m_allflow_dev_vec.pvnode_vec.size();++count)
	 	  		{
	 	  			j = m_allflow_dev_vec.pvnode_vec[count].i;
	 	  			B[j] = 0;
	 	  		}
			}
	 	  	n_u = 0;
	 	  	i_above = 0;
	 	    while(i_above <= i-1)
	 	  	{
	 	  	    count = 1;
	 	  	    while(count <= m_nusum_vec[i_above])
	 	        {
	 	  	        if(m_u_vec[n_u].j == i)//?
	 	  			{
	 	  			    Btemp = m_u_vec[n_u].value/m_d_vec[i_above];
						while(count <= m_nusum_vec[i_above])
						{
							j = m_u_vec[n_u].j;
							B[j] -= Btemp*m_u_vec[n_u].value;
							count++ ;
							n_u++;
						}
	 	  			    break;
	 	  			}
	 	  			else
	 	  			{
	 	  			    count++;
	 	  			    n_u++;
	 	  			}
				}
	 	  		
	 	  		i_above++;
	 	  	}
			float temp = B[i];
	 	  	Btemp = 1/B[i];
	 	  	m_d_vec[i] = Btemp;
	 	  	count = 0;
	 	  	for(j= i+1;j<m_node_num-1;j++)
	 	  	{
	 	  		    if(B[j] != 0)
	 	  		    {
						float  value= B[j];
	 	  				m_u_vec[n_u].value = B[j]*Btemp;
						temp = m_u_vec[n_u].value;
	 	  				m_u_vec[n_u].j = j;
	 	  				count++;
	 	  				n_u++;
	 	  		    }
			}
			m_nusum_vec[i] = count;
	 	}
	 }
	 for(i = 0; i<m_node_num;++i)
	 {
		 printf("i === %d,d === %f\n",i,m_d_vec[i]);
	 }
	 for(i = 0; i<m_node_num;++i)
	 {
		 printf("i === %d,u === %f\n",i,m_u_vec[i].value);
	 }
}

void CNetFlow::FormFactorTable1(int flag)
{
	int n_pv = 0;
	int n_u;//因子表上三角矩阵元素计数变量
	int count,i_above;//消去行号的计数变量
	int i_pv = m_allflow_dev_vec.pvnode_vec[0].i;
	int i,j;
	vector<float> B;
	float Btemp;
	B.resize(m_node_num);