📄 nflow_sub.cpp
字号:
for(I=0;I<2*lv_bus;I++)
{
B1[I]=0.0;
B2[I]=0.0;
}
if(MA==1)
{
//dP
//P
for(JSQ=0;JSQ<Ym[MX].lv_nz;JSQ++)
{//P的行,形成一行
//printf("Ym[%d].nz[%d].j==%d\n",MX,JSQ,Ym[MX].nz[JSQ].j);
J=2*Ym[MX].nz[JSQ].j;
B1[J]=matrix_jacb_newton[MX].EJU_CJPQ[JSQ*2];//Hij
B1[J+1]=matrix_jacb_newton[MX].EJU_CJPQ[JSQ*2+1];//Nij
B2[J]=matrix_jacb_newton[MX].EJL_CJPQ[JSQ*2];//Hji
B2[J+1]=matrix_jacb_newton[MX].EJL_CJPQ[JSQ*2+1];//Jji
}
B1[K+1]=matrix_jacb_newton[MX].DJF_CJPQ[0];//Nii,行,Hii是雅可比矩阵的真正对角元 并没有放入
B2[K+1]=matrix_jacb_newton[MX].DJF_CJPQ[1];//Jii,列
MA=0;
}
else if(MA==0)
{//Q行
// if(typ_bs[map_bus[MX].ibs]==BS_TYPE_PQ)//PQ 节点
if(map_bs[map_bus[MX].ibs].bs_type==BS_TYPE_PQ)//PQ 节点
{/* ! MX is a PQ bus*/
for(JSQ=0;JSQ<Ym[MX].lv_nz;JSQ++)
{
J=2*Ym[MX].nz[JSQ].j;
B1[J]=-matrix_jacb_newton[MX].EJU_CJPQ[JSQ*2+1];//Jij
B1[J+1]=matrix_jacb_newton[MX].EJU_CJPQ[JSQ*2];//Lij
B2[J]=-matrix_jacb_newton[MX].EJL_CJPQ[JSQ*2+1];//Nji
B2[J+1]=matrix_jacb_newton[MX].EJL_CJPQ[JSQ*2];//Lji
}
// MX=MX+1;
MA=1;
}
else
{
// MX=MX+1;
MA=1;
// printf("222222222222222\n");
continue;
}
}
}
if(flag==1)
{
for(i=0;i<lv_bus;i++)
{
// if(typ_bs[map_bus[MX].ibs]!=BS_TYPE_PQ)
if(map_bs[map_bus[i].ibs].bs_type!=BS_TYPE_PQ)
{/* ! BUS is a PV node*/
//PV节点，V部分置0,不参与计算
B1[2*i+1]=0.0;
B2[2*i+1]=0.0;
}
}
B1[2*(lv_bus-1)]=0.0;
B2[2*(lv_bus-1)]=0.0;
matrix_fact_newton[K].UD_FACT=matrix_jacb_newton[MX].DJ_CJPQ[MA];//P:Hii/Q:Lii,K行的对角元,MA作为标志，决定是选取Hii还是Lii
// printf("matrix_fact_newton[%d].UD_FACT==%f\n",K,matrix_fact_newton[K].UD_FACT);
//DEBUG
// printf("###################K==%d#########################################\n",K);
// printf("mx_jac[%d][%d]==%f\n",K,K,matrix_fact_newton[K].UD_FACT);
// for(JSQ=0;JSQ<lv_bus*2;JSQ++)
// {
// if(fabs(B1[JSQ])>0.00001)
// {
// printf("mx_jac[%d][%d]==%f \n",K,JSQ,B1[JSQ]);
// }
// }
//
// printf("\n");
// for(JSQ=0;JSQ<lv_bus*2;JSQ++)
// {
// if(fabs(B2[JSQ])>0.00001)
// {
// printf("mx_jac[%d][%d]==%f \n",JSQ,K,B2[JSQ]);
// }
// }
// printf("\n");
//DEBUG
for(I=0;I<K;I++)
{//
// printf("&&&&&&&&&&&&&&&&&&I=%d&&&&&&&&&&&&&&&&&&\n",I);
for(JSQ=0;JSQ<matrix_fact_newton[I].lv_nz;JSQ++)//NFD_FACT[]:行非对角元个数
{//I行的所有元
J=matrix_fact_newton[I].nz[JSQ].j-1;//非对角元列指
// printf("=================J==%d================\n",J);
if(J==K)
{
// printf("$$$$$$$$$$$$$$$$$$$$$$$$I=%d,lv_nz==%d,JSQ==%d$$$$$$$$$$$$$$$$$$$$$$$\n",I,matrix_fact_newton[I].lv_nz,JSQ);
// printf("=================JJ==%d================\n",J);
AI=matrix_fact_newton[I].nz[JSQ].UL_FACT;//下三角
BI=matrix_fact_newton[I].nz[JSQ].UU_FACT;//上三角
// printf("AI: UL_FACT[%d][%d]==%f, BI: UU_FACT[%d][%d]==%f\n",J,I,AI,I,J,BI);
//
// printf("UD_FACT[%d]==%f \n",K,matrix_fact_newton[K].UD_FACT);
matrix_fact_newton[K].UD_FACT=matrix_fact_newton[K].UD_FACT-AI*BI;//对角元
// printf("UD_FACT[%d]==%f \n\n",K,matrix_fact_newton[K].UD_FACT);
for(J1=JSQ+1;J1<matrix_fact_newton[I].lv_nz;J1++)
{
J=matrix_fact_newton[I].nz[J1].j-1;
// printf("#########J==%d#############\n",J);
// printf("B1[%d]=%f\n",J,B1[J]);
B1[J]=B1[J]-AI*matrix_fact_newton[I].nz[J1].UU_FACT;//行
//DEBUG
// printf("UU_FACT[%d][%d]==%f\n",I,J,matrix_fact_newton[I].nz[J1].UU_FACT);
// printf("B1[%d]=B1[%d]-AI*UU_FACT[%d][%d]\n",J,J,I,J);
// printf("B1[%d]==%f\n",J,B1[J]);
// printf("B2[%d]=%f\n",J,B2[J]);
//DEBUG
B2[J]=B2[J]-BI*matrix_fact_newton[I].nz[J1].UL_FACT;
//DEBUG
// printf("UL_FACT[%d][%d]==%f\n",J,I,matrix_fact_newton[I].nz[J1].UL_FACT);
// printf("B2[%d]=B2[%d]-BI*UL_FACT[%d][%d]\n",J,J,I,J);
// printf("B2[%d]==%f\n",J,B2[J]);
// printf("\n");
//DEBUG
}
break;
}
else if(J>K)
{
break;
}
}
}
// printf("#################################\n");
//DEBUG 07.07.09 matrix_fact_newton[K].lv_nz=0;
AI=matrix_fact_newton[K].UD_FACT;
if(fabs(AI)>0.00001)
{
AI=1.0/AI;
matrix_fact_newton[K].UD_FACT=AI;//对角元值
// printf("matrix_fact_newton[%d].UD_FACT==%f \n",K,matrix_fact_newton[K].UD_FACT);
for(J=K+1;J<2*lv_bus;J++)
{
if(fabs(B1[J])>0.00001||fabs(B2[J])>0.00001)
{
matrix_fact_newton[K].nz[matrix_fact_newton[K].lv_nz].UU_FACT=B1[J];
matrix_fact_newton[K].nz[matrix_fact_newton[K].lv_nz].UL_FACT=B2[J]*AI;//下三角元除对角元
matrix_fact_newton[K].nz[matrix_fact_newton[K].lv_nz].j=J+1;
// printf("j=%d, matrix_fact_newton[%d].nz[%d].UU_FACT=%f\n",matrix_fact_newton[K].nz[matrix_fact_newton[K].lv_nz].j-1,K,matrix_fact_newton[K].lv_nz,matrix_fact_newton[K].nz[matrix_fact_newton[K].lv_nz].UU_FACT);
// printf("j=%d, matrix_fact_newton[%d].nz[%d].UL_FACT=%f\n",matrix_fact_newton[K].nz[matrix_fact_newton[K].lv_nz].j-1,K,matrix_fact_newton[K].lv_nz,matrix_fact_newton[K].nz[matrix_fact_newton[K].lv_nz].UL_FACT);
matrix_fact_newton[K].lv_nz++;
}
}
}
}
}
else if(MA==1)
{
MA=0;
}
else
{
MA=1;
}
}
for(I=0;I<2*lv_bus;I++)
{
// printf("555:D[%d]==%f\n",I,matrix_fact_newton[I].UD_FACT);
for(J=0;J<matrix_fact_newton[I].lv_nz;J++)
{
matrix_fact_newton[I].nz[J].UU_FACT=matrix_fact_newton[I].nz[J].UU_FACT*matrix_fact_newton[I].UD_FACT;
}
}
// printf("\n\n");
// print_debug_info(7);
// print_debug_info(9);
// print_debug_info(21);
return RET_OK;
}
int qjj_newton(int NOWISLAN)
{
//按(dO1,dV1/V1,dO2,dV2/V2,...dOn,dVn/Vn)顺序排列变量
//三角分解法解方程
int I,J,JSQ;
// printf("OPT_FLAG=%d\n",OPT_FLAG);
// printf("111111111111111111\n");
// print_debug_info(10);
if(OPT_FLAG==0)
{
C[2*(lv_bus-1)]=0.0;//VO节点
for(I=lv_bus-NOV-1;I<lv_bus;I++)
{
// if(typ_bs[map_bus[I].ibs]!=BS_TYPE_PQ)//PV节点
if(map_bs[map_bus[I].ibs].bs_type!=BS_TYPE_PQ)//PV节点
{
C[2*I+1]=0.0;//dQ==0
}
}
// printf("11111111111111111\n\n");
//原始方程: YV=I
//导纳矩阵分解:Y=LDU,三角分解法解方程顺序: LZ=I,DX=Z,UV=X
//LZ=I方程
for(I=0;I<2*lv_bus;I++)
{
// printf("I==%d\n\n",I);
for(JSQ=0;JSQ<matrix_fact_newton[I].lv_nz;JSQ++)
{
J=matrix_fact_newton[I].nz[JSQ].j-1;
// {
// printf("J==%d\n",J);
// printf("C[%d]=%f,C[%d]==%f\n",I,C[I],J,C[J]);
// printf("L[%d][%d]==%f\n",J,I,matrix_fact_newton[I].nz[JSQ].UL_FACT);
// printf("C[%d]=C[%d]-C[%d]*L[%d][%d]\n\n",J,J,I,J,I);
// }
C[J]=C[J]-C[I]*matrix_fact_newton[I].nz[JSQ].UL_FACT;//按列
// printf("111:C[%d]=%f\n",J,C[J]);
}
}
// printf("22222222222222222222\n");
// print_debug_info(10);
// printf("22222222222222222222\n");
for(I=2*lv_bus-1;I>-1;I--)
{
//DX=Z
// printf("D[%d]=%f\n",I,matrix_fact_newton[I].UD_FACT);
// printf("C[%d]=%f,D%d=%f\n",I,C[I],I,matrix_fact_newton[I].UD_FACT);
C[I]=C[I]*matrix_fact_newton[I].UD_FACT;
// printf("222:C[%d]=%f\n",I,C[I]);
//UV=X
for(JSQ=0;JSQ<matrix_fact_newton[I].lv_nz;JSQ++)
{
J=matrix_fact_newton[I].nz[JSQ].j-1;
//if(I==484&&IT==1)
// {
// printf("111:C[%d]=%f,C[%d]==%f\n",I,C[I],J,C[J]);
// printf("U[%d][%d]==%f\n",I,J,matrix_fact_newton[I].nz[JSQ].UU_FACT);
// printf("C[%d]=C[%d]-C[%d]*U[%d][%d]\n",I,I,J,I,J);
// }
C[I]=C[I]-C[J]*matrix_fact_newton[I].nz[JSQ].UU_FACT;
//if(I==484&&IT==1)
// printf("333:C[%d]=%f\n\n",I,C[I]);
}
}
// printf("4444444444444444\n\n");
C[2*(lv_bus-1)]=0.0;
// for(I=lv_bus-NOV-1;I<lv_bus;I++)
for(I=0;I<lv_bus;I++)
{
// if(typ_bs[map_bus[I].ibs]!=BS_TYPE_PQ)
if(map_bs[map_bus[I].ibs].bs_type!=BS_TYPE_PQ)
{
C[2*I+1]=0.0;
// printf("C[%d]==%f\n",2*I+1,C[2*I+1]);
}
}
// printf("5555555555555555555555\n\n");
//计算新的电压/角度值
for(I=0;I<lv_bus;I++)
{
map_bus[I].v_pu=map_bus[I].v_pu-C[2*I+1]*map_bus[I].v_pu;
map_bus[I].a_pu=map_bus[I].a_pu-C[2*I];
}
}
// else if(OPT_FLAG==1)
// {//NO DEBUG
// for(I=0;I<2*lv_bus;I++)
// {
// for(JSQ=0;JSQ<matrix_fact_newton[I].lv_nz;JSQ++)
// {
// J=matrix_fact_newton[I].nz[JSQ].j-1;
// B[J]=B[J]-B[I]*matrix_fact_newton[I].nz[JSQ].UU_FACT;
//
// }
// }
// for(I=2*lv_bus-1;I>-1;I--)
// {
// B[I]=B[I]*matrix_fact_newton[I].UD_FACT;
// for(JSQ=0;JSQ<matrix_fact_newton[I].lv_nz;JSQ++)
// {
//
// J=matrix_fact_newton[I].nz[JSQ].j;
// B[I]=B[I]-B[J-1]*matrix_fact_newton[I].nz[JSQ].UL_FACT;
// }
// }
// }
// printf("33333333333333333333333333\n");
print_debug_info(10);
return RET_OK;
}
int pf_bus_to_bs(int NOWISLAND)
{
int I;
for(I=0;I<lv_bs;I++)
{
if(map_bs[I].AC_island==NOWISLAND)
{
map_bs[I].v_ym=map_bus[map_bs[I].ibus].v_pu*map_vt[map_vl[map_bs[I].vl_no].vt_no].vbase;
map_bs[I].a_ym=map_bus[map_bs[I].ibus].a_pu*180.0/3.14159;
}
}
pf_bus(NOWISLAND);
pf_branch(NOWISLAND);
pf_inj(NOWISLAND);
// for(I=0;I<lv_map_AClnbrh;I++)
// {
// if(map_AClnbrh[I].iln2>-1&&map_bs[map_nd[map_AClnbrh[I].ind].ibs-1].AC_island==NOWISLAND)
// {
// if(fabs(map_ACln2[map_AClnbrh[I].iln2].Pi-map_ACln2[map_AClnbrh[I].iln2].p/100.0)>0.1)
// {
// printf("ln_name_i==%s, Pij==%f, dlt_p==%f\n",map_ACln[map_AClnbrh[I].ACln_no].name,map_ACln2[map_AClnbrh[I].iln2].Pi,map_ACln2[map_AClnbrh[I].iln2].Pi-map_ACln2[map_AClnbrh[I].iln2].p/100.0);
// }
//
// if(fabs(map_ACln2[map_AClnbrh[I].iln2].Qi-map_ACln2[map_AClnbrh[I].iln2].q/100.0)>0.1)
// {
// printf("ln_name_i==%s, Qij==%f, dlt_q==%f\n",map_ACln[map_AClnbrh[I].ACln_no].name,map_ACln2[map_AClnbrh[I].iln2].Qi,map_ACln2[map_AClnbrh[I].iln2].Qi-map_ACln2[map_AClnbrh[I].iln2].q/100.0);
// }
// }
//
// if(map_AClnbrh[I].zln2>-1&&map_bs[map_nd[map_AClnbrh[I].znd].ibs-1].AC_island==NOWISLAND)
// {
// if(fabs(map_ACln2[map_AClnbrh[I].zln2].Pi-map_ACln2[map_AClnbrh[I].zln2].p/100.0)>0.1)
// {
// printf("ln_name_z==%s, Pij==%f, dlt_p==%f\n",map_ACln[map_AClnbrh[I].ACln_no].name,map_ACln2[map_AClnbrh[I].zln2].Pi,map_ACln2[map_AClnbrh[I].zln2].Pi-map_ACln2[map_AClnbrh[I].zln2].p/100.0);
// }
//
// if(fabs(map_ACln2[map_AClnbrh[I].iln2].Qi-map_ACln2[map_AClnbrh[I].iln2].q/100.0)>0.1)
// {
// printf("ln_name_z==%s, Qij==%f, dlt_q==%f\n",map_ACln[map_AClnbrh[I].ACln_no].name,map_ACln2[map_AClnbrh[I].zln2].Qi,map_ACln2[map_AClnbrh[I].zln2].Qi-map_ACln2[map_AClnbrh[I].zln2].q/100.0);
// }
// }
// }
// for(I=0;I<lv_map_ACxfbrh;I++)
// {
// if(map_ACxfbrh[I].ixf>-1&&map_bs[map_nd[map_ACxfbrh[I].ind].ibs-1].AC_island==NOWISLAND)
// {
// if(fabs(map_ACxf[map_ACxfbrh[I].ixf].Pi-map_ACxf[map_ACxfbrh[I].ixf].p/100.0)>0.1)
// {
// printf("xf_name==%s, st_name==%s, Pij==%f, dlt_p==%f\n",map_ACxf[map_ACxfbrh[I].ixf].name,map_st[map_ACxf[map_ACxfbrh[I].ixf].st_no].name,map_ACxf[map_ACxfbrh[I].ixf].Pi,map_ACxf[map_ACxfbrh[I].ixf].Pi-map_ACxf[map_ACxfbrh[I].ixf].p/100.0);
// }
//
// if(fabs(map_ACln2[map_AClnbrh[I].iln2].Qi-map_ACln2[map_AClnbrh[I].iln2].q/100.0)>0.1)
// {
// printf("xf_name==%s, st_name==%s, Qij==%f, dlt_q==%f\n",map_ACxf[map_ACxfbrh[I].ixf].name,map_st[map_ACxf[map_ACxfbrh[I].ixf].st_no].name,map_ACxf[map_ACxfbrh[I].ixf].Qi,map_ACxf[map_ACxfbrh[I].ixf].Qi-map_ACxf[map_ACxfbrh[I].ixf].q/100.0);
// }
//
//
// }
// }
return RET_OK;
}
int dy_xi_dpx(int NOWISLAND)
{
float MS;
int I,IK;
//C[I]的顺序：I,DPi,I+1,DQi
if(IT<1)
{//首次计算
MS=0.0;
for(I=0;I<lv_bus;I++)
{
map_bus[I].dyn_xi=0.0;
}
// if(SYS_DIM==1)
// {
// for(I=0;I<lv_map_ACun;I++)
// {
// if(map_bs[map_nd[map_ACun[I].ind].ibs-1].AC_island==NOWISLAND)
// {
// if(netflt_->XD_UN[i]>0.0001)
// {//一般不用这种方式计算系数
// IK=netpf_->JA_BS[netbus_->I$BS_UN[i]-1]-1;
// pfpart1_->XI_BUS[IK]=pfpart1_->XI_BUS[IK]+1.0/netflt_->XD_UN[i];
// MS=MS+1.0/netflt_->XD_UN[i];//XD_UN??
// }
// }
// }
// }
// else
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -