📄 nlf.cpp
字号:
int i,sum=0,count;
gene_yymatrix[1].addres=1;
for(i=1;i<=N;i++)
{
if(i==NS)
{
gene_yymatrix[2*i].addres=gene_yymatrix[2*i-1].addres;
gene_yymatrix[2*i+1].addres=gene_yymatrix[2*i].addres;
continue;
}
int n,j,end,list;
for(j=1;j<i;j++)//先求取某一节点p，q或v2对f和e的偏导形成的工作数组//考虑1---i节点
{
if(j==NS)continue;//排除平衡节点
for(end=yy_matrix[j+1].addres-1,n=yy_matrix[j].addres;n<=end;n++)
{
if(y_matrix[n].number==i)
{
UP[2*j-1]+=y_matrix[n].Bij*voltage[i].e-y_matrix[n].Gij*voltage[i].f;//f排在前面，e排在后面
UP[2*j]+=-(y_matrix[n].Gij*voltage[i].e+y_matrix[n].Bij*voltage[i].f);
if(pv[i]!=1)
{
UQ[2*j-1]+=-UP[2*j];
UQ[2*j]+=UP[2*j-1];
}
// continue;考虑双回路
}
}
} //考虑自身倒数
UP[2*i-1]=-ab[i].bii+(yy_matrix[i].Bii*voltage[i].e-yy_matrix[i].Gii*voltage[i].f);//p对f的倒数
UP[2*i]=-ab[i].aii-(yy_matrix[i].Gii*voltage[i].e+yy_matrix[i].Bii*voltage[i].f);//p对e的倒数
if(pv[i]!=1)
{
UQ[2*i-1]=-ab[i].aii+(yy_matrix[i].Gii*voltage[i].e+yy_matrix[i].Bii*voltage[i].f);
UQ[2*i]=ab[i].bii+(yy_matrix[i].Bii*voltage[i].e-yy_matrix[i].Gii*voltage[i].f);
}
for(end=yy_matrix[i+1].addres-1,n=yy_matrix[i].addres;n<=end;n++)//考虑i--N的倒数
{
list=y_matrix[n].number;
if(list==NS)continue;//排除平衡节点
UP[2*list-1]+=y_matrix[n].Bij*voltage[i].e-y_matrix[n].Gij*voltage[i].f;
UP[2*list]+=-(y_matrix[n].Gij*voltage[i].e+y_matrix[n].Bij*voltage[i].f);
if(pv[i]!=1)
{
UQ[2*list-1]+=-UP[2*list];
UQ[2*list]+=UP[2*list-1];//考虑双回路
}
}
if(pv[i]==1)//求pv节点的v2对f和e的偏导
{
UQ[2*i-1]=-2*voltage[i].f;
UQ[2*i]=-2*voltage[i].e;
}
//再消去已生成的两行工作数组，再稀疏存储形成因子表，同时对dpqv右端向量进行前代处理
for(j=1;j<=2*i-2;j++)//对p对f和e的偏导数消元
{
if(int((j+1)/2)==NS)continue;
if(UP[j]!=0)
{
for(end=gene_yymatrix[j+1].addres-1,n=gene_yymatrix[j].addres;n<=end;n++)
{
list=gene_ymatrix[n].number;
UP[list]-=UP[j]*gene_ymatrix[n].data;//消去工作数组
}
if(j%2==1)
dpqv[i].dp-=UP[j]*dpqv[int(j+1)/2].dp;//消去右端向量
else
if(pv[int(j+1)/2]==1)
dpqv[i].dp-=UP[j]*dpqv[int(j+1)/2].dv2;
else
dpqv[i].dp-=UP[j]*dpqv[int(j+1)/2].dq;
}
}
//此行稀疏存储p对f和e的偏导数消元后的工作数组
count=0;
for(j=2*i;j<=2*N;j++)
{
if(int(j+1)/2==NS)continue;
if(fabs(UP[j])>0)
{
count++;
sum++;
gene_ymatrix[sum].data=UP[j]/UP[2*i-1];//标准化
gene_ymatrix[sum].number=j;
}
}
gene_yymatrix[2*i-1].data=UP[2*i-1];
gene_yymatrix[2*i].addres=gene_yymatrix[2*i-1].addres+count;
dpqv[i].dp=dpqv[i].dp/UP[2*i-1];//标准化右端向量
for(j=1;j<=2*i-1;j++)//对q或v2对f和e的偏导数消元
{
if(int(j+1)/2==NS)continue;
if(UQ[j]!=0)
{
for(end=gene_yymatrix[j+1].addres-1,n=gene_yymatrix[j].addres;n<=end;n++)
{
list=gene_ymatrix[n].number;
UQ[list]-=gene_ymatrix[n].data*UQ[j];
}
if(pv[i]==1)
{
if(j%2==1)
dpqv[i].dv2-=UQ[j]*dpqv[int(j+1)/2].dp;//消去右端向量
else
if(pv[int(j+1)/2]==1)
dpqv[i].dv2-=UQ[j]*dpqv[int(j+1)/2].dv2;
else
dpqv[i].dv2-=UQ[j]*dpqv[int(j+1)/2].dq;
}
else
{
if(j%2==1)
dpqv[i].dq-=UQ[j]*dpqv[int(j+1)/2].dp;//消去右端向量
else
if(pv[int(j+1)/2]==1)
dpqv[i].dq-=UQ[j]*dpqv[int(j+1)/2].dv2;
else
dpqv[i].dq-=UQ[j]*dpqv[int(j+1)/2].dq;
}
}
}
count=0;
for(j=2*i+1;j<=2*N;j++)//稀疏存储dq或dv2对f和e偏导数消元后的工作数组
{
if(int(j+1)/2==NS)continue;
if(fabs(UQ[j])>0)
{
count++;
sum++;
gene_ymatrix[sum].data=UQ[j]/UQ[2*i];//标准化
gene_ymatrix[sum].number=j;
}
}
gene_yymatrix[2*i].data=UQ[2*i];
gene_yymatrix[2*i+1].addres=gene_yymatrix[2*i].addres+count;
if(pv[i]==1)
dpqv[i].dv2=dpqv[i].dv2/UQ[2*i];//标准化右端向量
else
dpqv[i].dq=dpqv[i].dq/UQ[2*i];
for(j=1;j<=2*N;j++)//工作数组清零
{
UP[j]=0;
UQ[j]=0;
}
}
free(UP);
free(UQ);
}
void nlf::back_solve()
{
int j,n,end,list;
for(j=N;j>=1;j--)
{
if(j==NS)continue;
for(end=gene_yymatrix[2*j+1].addres-1,n=gene_yymatrix[2*j].addres;n<=end;n++)
{
if(pv[j]==1)
{
list=gene_ymatrix[n].number;
if(list%2==1)
dpqv[j].dv2-=gene_ymatrix[n].data*dpqv[int(list+1)/2].dp;
else
if(pv[int(list+1)/2]==1)
dpqv[j].dv2-=gene_ymatrix[n].data*dpqv[int(list+1)/2].dv2;
else
dpqv[j].dv2-=gene_ymatrix[n].data*dpqv[int(list+1)/2].dq;
}
else
{
list=gene_ymatrix[n].number;
if(list%2==1)
dpqv[j].dq-=gene_ymatrix[n].data*dpqv[int(list+1)/2].dp;
else
if(pv[int(list+1)/2]==1)
dpqv[j].dq-=gene_ymatrix[n].data*dpqv[int(list+1)/2].dv2;
else
dpqv[j].dq-=gene_ymatrix[n].data*dpqv[int(list+1)/2].dq;
}
}
for(end=gene_yymatrix[2*j].addres-1,n=gene_yymatrix[2*j-1].addres;n<=end;n++)
{
list=gene_ymatrix[n].number;
if(list%2==1)
dpqv[j].dp-=gene_ymatrix[n].data*dpqv[int(list+1)/2].dp;
else
if(pv[int(list+1)/2]==1)
dpqv[j].dp-=gene_ymatrix[n].data*dpqv[int(list+1)/2].dv2;
else
dpqv[j].dp-=gene_ymatrix[n].data*dpqv[int(list+1)/2].dq;
}
}
}
void nlf::amend_voltage()
{
int i;
for(i=1;i<=N;i++)
{
if(i==NS)continue;
voltage[i].f-=dpqv[i].dp;
if(pv[i]==1)
voltage[i].e-=dpqv[i].dv2;
else
voltage[i].e-=dpqv[i].dq;
}
}
void nlf::BeginPowerflow()
{
fscanf(in,"%d %d %d %d",&N,&NS,&number_branch,&number_ground);
branch=(BRANCH *)calloc((number_branch+1), sizeof(BRANCH));//原始数据存储申请
if(!branch)EXIT;
genrt=(GENRT *)calloc((N+1),sizeof(GENRT));
if(!genrt)EXIT;
load=(LOAD *)calloc((N+1),sizeof(LOAD));
if(!load)EXIT;
ground=(GROUND *)calloc((number_ground+1),sizeof(GROUND));
if(!ground)EXIT;
ppv=(int *)calloc((N+1),sizeof(int));
if(!ppv)EXIT;
node=(NODE *)calloc((N+1),sizeof(NODE));
if(!node)EXIT;
pv=(int *)calloc((N+1),sizeof(int)); //节点优化数据
if(!pv)EXIT;
o_n=(int *)calloc((N+1),sizeof(int));
if(!o_n)EXIT;
n_o=(int *)calloc((N+1),sizeof(int));
if(!n_o)EXIT;
y_matrix=(Y_MATRIX *)calloc((number_branch+1),sizeof(Y_MATRIX));//申请导纳阵存储空间
if(!y_matrix)EXIT;
yy_matrix=(YY_MATRIX *)calloc((N+2),sizeof(YY_MATRIX));
if(!yy_matrix)EXIT;
gene_ymatrix=(GENE_YMATRIX *)calloc((20*number_branch+1),sizeof(GENE_YMATRIX));//申请修正矩阵存储空间
if(!gene_ymatrix)EXIT;
gene_yymatrix=(GENE_YYMATRIX *)calloc((2*N+2),sizeof(GENE_YYMATRIX));
if(!gene_yymatrix)EXIT;
dpqv=(PQV *)calloc((N+1),sizeof(PQV));//申请右端向量存储空间
if(!dpqv)EXIT;
voltage=(VOLTAGE *)calloc((N+1),sizeof(VOLTAGE));//申请电压存储空间
if(!voltage)EXIT;
ab=(AB*)calloc((N+1),sizeof(AB));//申请修正含量aii和bii所用存储空间
if(!ab)EXIT;
readdata();
node_order();
re_order();
branch_rank();
form_ymatrix();
initia_pqv();
}
void nlf::iterative()
{
bool mark=true;
int j;
double big1,big2;
iter_number=0;
while(mark==true)
{
initia_dpqv();
gene_matrix();
back_solve();
amend_voltage();
iter_number++;
if(pv[1]==1)
big1=(dpqv[1].dp>dpqv[1].dv2)?dpqv[1].dp:dpqv[1].dv2;
else
big1=(dpqv[1].dp>dpqv[1].dq)?dpqv[1].dp:dpqv[1].dq;
for(j=2;j<=N;j++)
{
if(pv[j]==1)
big2=(dpqv[j].dp>dpqv[j].dv2)?dpqv[j].dp:dpqv[j].dv2;
else
big2=(dpqv[j].dp>dpqv[j].dq)?dpqv[j].dp:dpqv[j].dq;
if(big1<big2)
big1=big2;
}
fprintf(out,"%-12s%-6d%-12s%-15.10lf\n","迭代次数:",iter_number,"最大修正值:",big1);
if(big1>EPS&&iter_number<8)
mark=true;
else
mark=false;
}
}
void nlf::result()
{
int i,new_node;
PointData module, angle,loss=0.0;
fprintf(out,"\n%-12s%-6d%-12s%-15.10lf\n\n","迭代次数:",iter_number,"迭代精度:",EPS);
fprintf(out,"%-12s%-20s%-21s%-22s%-20s\n\n","节点号","实部电压"," 虚部电压","幅值","角度");
for(i=1;i<=N;i++)
{
new_node=o_n[i];
module=sqrt(voltage[new_node].e*voltage[new_node].e+voltage[new_node].f*voltage[new_node].f);
angle=atan(voltage[new_node].f/voltage[new_node].e)/3.1415926*180;
fprintf(out,"%-10d%-20.6lf%-20.6lf%-20.6lf%-20.6lf\n",i,voltage[new_node].e,
voltage[new_node].f,module,angle);
}
fprintf(out,"\n\n%-10s%-10s%-15s%-15s%-15s%-15s%-15s%-15s\n\n", "i","j","Pij","Qij","Pji","Qji","lossp","lossq");
for(i=1;i<=number_branch;i++)
{
int inode,jnode;
PointData Pij,Qij,Pji,Qji,LPij,LQij;
if(branch[i].bk<0)
{
inode=branch[i].i_node;
jnode=branch[i].j_node;
if(inode<0)
{
inode=-inode;
voltage[inode].e=-voltage[inode].e/branch[i].bk;
voltage[inode].f=-voltage[inode].f/branch[i].bk;
}
if(jnode<0)
{
jnode=-jnode;
voltage[jnode].e=-voltage[jnode].e/branch[i].bk;
voltage[jnode].f=-voltage[jnode].f/branch[i].bk;
}
Pij=voltage[inode].e*(y_matrix[i].Gij*voltage[jnode].e-y_matrix[i].Bij*voltage[jnode].f)
+voltage[inode].f*(y_matrix[i].Gij*voltage[jnode].f+y_matrix[i].Bij*voltage[jnode].e)
+(voltage[inode].e*voltage[inode].e+voltage[inode].f*voltage[inode].f)*
(-y_matrix[i].Gij);
Qij=voltage[inode].f*(y_matrix[i].Gij*voltage[jnode].e-y_matrix[i].Bij*voltage[jnode].f)
-voltage[inode].e*(y_matrix[i].Gij*voltage[jnode].f+y_matrix[i].Bij*voltage[jnode].e)
-(voltage[inode].e*voltage[inode].e+voltage[inode].f*voltage[inode].f)*
(-y_matrix[i].Bij);
inode=abs(branch[i].j_node);//交换节点
jnode=abs(branch[i].i_node);
Pji=voltage[inode].e*(y_matrix[i].Gij*voltage[jnode].e-y_matrix[i].Bij*voltage[jnode].f)
+voltage[inode].f*(y_matrix[i].Gij*voltage[jnode].f+y_matrix[i].Bij*voltage[jnode].e)
+(voltage[inode].e*voltage[inode].e+voltage[inode].f*voltage[inode].f)*
(-y_matrix[i].Gij);
Qji=voltage[inode].f*(y_matrix[i].Gij*voltage[jnode].e-y_matrix[i].Bij*voltage[jnode].f)
-voltage[inode].e*(y_matrix[i].Gij*voltage[jnode].f+y_matrix[i].Bij*voltage[jnode].e)
-(voltage[inode].e*voltage[inode].e+voltage[inode].f*voltage[inode].f)*
(-y_matrix[i].Bij);
LPij=Pij+Pji;
LQij=Qij+Qji;
inode=branch[i].i_node;
jnode=branch[i].j_node;
if(inode<0)
{
inode=-inode;
voltage[inode].e=-voltage[inode].e*branch[i].bk;
voltage[inode].f=-voltage[inode].f*branch[i].bk;
}
if(jnode<0)
{
jnode=-jnode;
voltage[jnode].e=-voltage[jnode].e*branch[i].bk;
voltage[jnode].f=-voltage[jnode].f*branch[i].bk;
}
}
else
{
inode=branch[i].i_node;
jnode=branch[i].j_node;
Pij=voltage[inode].e*(y_matrix[i].Gij*voltage[jnode].e-y_matrix[i].Bij*voltage[jnode].f)
+voltage[inode].f*(y_matrix[i].Gij*voltage[jnode].f+y_matrix[i].Bij*voltage[jnode].e)
+(voltage[inode].e*voltage[inode].e+voltage[inode].f*voltage[inode].f)*
(-y_matrix[i].Gij);
Qij=voltage[inode].f*(y_matrix[i].Gij*voltage[jnode].e-y_matrix[i].Bij*voltage[jnode].f)
-voltage[inode].e*(y_matrix[i].Gij*voltage[jnode].f+y_matrix[i].Bij*voltage[jnode].e)
-(voltage[inode].e*voltage[inode].e+voltage[inode].f*voltage[inode].f)*
(-y_matrix[i].Bij+branch[i].bk);
inode=branch[i].j_node;//交换节点
jnode=branch[i].i_node;
Pji=voltage[inode].e*(y_matrix[i].Gij*voltage[jnode].e-y_matrix[i].Bij*voltage[jnode].f)
+voltage[inode].f*(y_matrix[i].Gij*voltage[jnode].f+y_matrix[i].Bij*voltage[jnode].e)
+(voltage[inode].e*voltage[inode].e+voltage[inode].f*voltage[inode].f)*
(-y_matrix[i].Gij);
Qji=voltage[inode].f*(y_matrix[i].Gij*voltage[jnode].e-y_matrix[i].Bij*voltage[jnode].f)
-voltage[inode].e*(y_matrix[i].Gij*voltage[jnode].f+y_matrix[i].Bij*voltage[jnode].e)
-(voltage[inode].e*voltage[inode].e+voltage[inode].f*voltage[inode].f)*
(-y_matrix[i].Bij+branch[i].bk);
LPij=Pij+Pji;
LQij=Qij+Qji;
}
inode=abs(branch[i].i_node);
jnode=abs(branch[i].j_node);
fprintf(out,"%-10d%-10d%-15f%-15f%-15f%-15f%-15f%-15f\n", n_o[inode],n_o[jnode],Pij,Qij,Pji,Qji,LPij,LQij);
loss+=LPij;
}
fprintf(out,"\n\n%10s%f\n","全网有功损耗",loss);
}
void nlf::free_point()
{
free(pv);
free(ppv);
free(o_n);
free(n_o);
free(branch);
free(genrt);
free(load);
free(ground);
free(voltage);
free(dpqv);
free(y_matrix);
free(yy_matrix);
free(gene_ymatrix);
free(gene_yymatrix);
free(ab);
}
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -