【正文】 newsort[i]){idk[j] = i。} } } void initial(){//**** 本函數(shù)進(jìn)行初始化工作 ****//int i, k1。//輸入原始數(shù)據(jù)oldtonew()。for(i = 1。nnew[i] = i。//節(jié)點(diǎn)優(yōu)化mpj = nnew[mpj]。// sort the r,x and b for(i = 1。idk[i] = nnew[k1]。i {v[i] = v00。} // (v00)，電壓相角初值都為0。i yds[i] = 0。i {k1 = nnew[nob[i]]。} for(i = 1。newval(pg)。newval(pl)。newval(v0)。i v[i] = v0[i]。//nobt[] = 0: pq節(jié)點(diǎn)，v0[]存放的是最后一個節(jié)點(diǎn)數(shù)據(jù)，//對于pq節(jié)點(diǎn)，即為系統(tǒng)平均電壓值。//輸出af[]、v00和節(jié)點(diǎn)排序后的支路、節(jié)點(diǎn)和//接地電抗數(shù)據(jù)(僅僅查看中間結(jié)果)ya0()。(P407)}void printo(){ //**** 輸出af[]、v00和節(jié)點(diǎn)排序后的支路、節(jié)點(diǎn)和接地電抗數(shù)據(jù) ****//int i。fprintf(fp2, “n %%%”, af[0], af[1], v00)。39。fprintf(fp2, “nn*******ZLB*******n”)。i {fprintf(fp2, “n”)。fprintf(fp2, “%%%”, zr[i], zx[i], zyk[i])。39。fprintf(fp2, “nn*******BUS*******n”)。i {fprintf(fp2, “n”)。fprintf(fp2, “%%%%%”, pg[i], qg[i], pl[i], ql[i], v0[i])。39。fprintf(fp2,“nn******DKK******n”)。i {fprintf(fp2, “n”)。} } void dataio(){//**** 系統(tǒng)數(shù)據(jù)初始化 ****//int i。af[1] = 。v00 = 。scanf(“%s”, inname)。printf(“nplease output the name of data filen”)。fp2 = fopen(outname, “w”)。n, amp。nb, amp。// the number of node ,branches, node fscanf(fp1, “%lf %d %d %d %d”, amp。kk2, amp。bnsopton, amp。//precision, swing node,sort the node,iteration numbers for(i = 1。izl[i], amp。fscanf(fp1, “%lf %lf %lf ”, amp。zx[i], amp。} for(i = 1。nob[i], amp。fscanf(fp1, “%lf %lf %lf %lf %lf”, amp。qg[i], amp。ql[i], amp。} for(i = 1。idk[i], amp。} fclose(fp1)。int i, t。yy1()。//form the B matrix of P0 iteration y2()。//form the B matrix of QV iteration t = 0。kg[0] = kg[1] = 1。// calculating the power bbhl(kq0)。elseprintf(“Q: %dt%dt%fn”, t, dsd, dsm)。if(kq0 == 0)calc(iu1, u1, di1, nfd1, b)。for(i = 1。} } else kg[kq0] = 0。else{kq0 = 0。} if(t it1)break。fprintf(fp2, “n%s%d”, “times = ”, t)。// recover the data if sorted // newtoold()。// node data branch_output()。39。printc(39。, 78)。}void newval(double* aa){//**** 本函數(shù)將舊號換成新號 ****//int i, k1。i b[i] = 。i {k1 = nnew[nob[i]]。} for(i = 1。}void yzb(int t, int* iu, double* u, double* di, int* nfd){ //**** 本函數(shù)求因子表 ****////參數(shù)1為標(biāo)志(t=0 求B39。39。int jj, jj1, jj2, im, x, fd[NS]。nfd[1] = 1。i //nobt[] 存放的是節(jié)點(diǎn)類型, 0: pq節(jié)點(diǎn),1: pv節(jié)點(diǎn)。amp。j //|b[j] = 。//|if((kk2 == 0)amp。(t == 1)amp。(nobt[i]!=1))// 存在(t == 1)的情況，不多余。//af[1] i1 = ydz[i]。for(im = 1。jj2 = nfd[im + 1]ai * u[k]。}}}x = nfd[i]。ai = di[i]。i1 = i + 1。j {if(fabs(b[j]) ){u[x] = b[j] * ai。k++。}}fd[i] = k。di[i] = 。} fprintf(fp2, “n********U*********”)。i if(i % 3 == 1)fprintf(fp2, “n”)。} fprintf(fp2, “n********DI********”)。}void printf1(double* aa, int n){ //**** 本函數(shù)輸出aa[i],i=1n ****// int i。i {if(i % 5 == 1)fprintf(fp2, “n”)。} fprintf(fp2, “nn”)。(P417圖F19)double bi。for(i = 1。i1 = nfd[i]。for(j = i1。j ++){k = iu[j]。for(j = i2。j){k = iu[j]。i ++){ip = i。k3 = iabs(jzl[i])。j {k2 = iabs(izl[j])。if(k2 k1 ||(k2 == k1 amp。 k4 k3)){ip = j。k3 = k4。izl[i], amp。}} }void bnsopt(){ iswap(amp。jzl[ip])。zr[i], amp。swap(amp。zx[ip])。zyk[i], amp。//**** 節(jié)點(diǎn)優(yōu)化 ****//int ii1, ii2, zls2, nomax。int temp。for(i = 1。//先清零。for(i = 1。old[i + zls] = iabs(jzl[i])。for(i = 1。for(j = k1。old[i], amp。//交換整數(shù)old[i]、old[j]。}nomax = old[zls2]。 l = 1。i {ii1 = old[l]。for(j = l。if(ii1!= ii2){l = j。}nnew[i] ++。i //即為第 i 行的非零元素個數(shù)。k if(kq0!= 0){ai =yb[k]。}//yb[]存放導(dǎo)納陣非對角元的虛部(yg[] + j*yb[])else{ai = yg[k]。}//yg[]存放導(dǎo)納陣非對角元的實(shí)部(yg[] + j*yb[])j = iy[k]。vj = vi * v[j]。ci = va[i]bi。(程序框圖見P423)****//int i。double pl0, pg0。// dsm 即為最大功率誤差。i vi = v[i]。vj = v0[i]。v0[]存放的是最后一個節(jié)點(diǎn)數(shù)據(jù)。// nopt[] = 0: pq節(jié)點(diǎn)，對于pq節(jié)點(diǎn)，即為系統(tǒng)平均電壓值。if(kq0 == 0){ pl0 = pl[i]。}// 負(fù)荷p,發(fā)電機(jī)pelse{pl0 = ql[i]。}// 負(fù)荷q,發(fā)電機(jī)q if(kk2 == 0)pl0 = pl0 *((vipl0。amp。//pv節(jié)點(diǎn)(nobt[] =1)和平衡節(jié)點(diǎn)(mpj)不參與 //求最大功率誤差和常數(shù)項(xiàng)的運(yùn)算if(((kq0!= 1)||(nobt[i]!=1))amp。(i!= mpj)){ if(fabs(b[i]) fabs(dsm)){dsm = b[i]。}// dsm 即為最大功率誤差,dsd存放其對應(yīng)的節(jié)點(diǎn)號。// 計算修正方程式的常數(shù)項(xiàng)。// pv節(jié)點(diǎn)(nobt[] =1)和平衡節(jié)點(diǎn)(mpj)不參與// 求最大功率誤差和常數(shù)項(xiàng)的運(yùn)算。int i, j, oldnumber。+39。fprintf(fp2, “n%5s%8s%10s%11s%11s%11s%11sn”, “I”，“V”, “CA”, “PL”, “QL”, “PG”, “QG”)。dsd = 1。i j = nnew[i]。//轉(zhuǎn)化為相應(yīng)舊號ci = va[j] * 。vi = v[j]。dsd = j。fprintf(fp2, “n%5d%%”, oldnumber, vi, ci)。} printc(39。, 72)。(程序框圖見p428 F117)int ii, jj。int i, j, l。****// dph = 。dqh = 。fprintf(fp2, “n%5s%5s%10s%12s%12s%12sn ”, “I”, “J”, “PIJ”，“QIJ”, “PJI”, “QJI”)。i {j = idk[i]。} for(l = 1。// izl[]: 支路左節(jié)點(diǎn)號。// jzl[]: 支路右節(jié)點(diǎn)號。j=nnew[jj]。ii = newsort[ii]。//轉(zhuǎn)化為相應(yīng)舊號 r = zr[l]。yk = zyk[l]。//v[]: 電壓幅值。//va[]: 電壓相角。cj = va[j]。fi = vi * sin(ci)。ej = vj * cos(cj)。// ej: 支路右節(jié)點(diǎn)電壓實(shí)部，fj: 支路右節(jié)點(diǎn)電壓虛部。{ if(izl[l] 0){ei = ei / yk。}// yk=zyk[l]else{ej = ej / yk。}yk = 。zf = r * r + x * x。fir =(df * rfir * ei。qji = qji1。k1 = iabs(izl[i])。for(j = i + 1。k4 = iabs(jzl[j])。amp。k1 = k2。}}if(i!= ip){iswap(amp。izl[ip])。jzl[i], amp。swap(amp。zr[ip])。zx[i], amp。swap(amp。zyk[ip])。i k = nob[i]。printf(“nob[%d] = %dn”, i, nob[i])。i for(j = i+1。nob[i] = nob[j]。}} } printf(“nob[%d] = %dn”, i, nob[i])。double r, x, yk, zf, gij, bij。for(i = 1。bii[i] = 。for(i = 1。yb[i] = 。ll = 1。l i1 = izl[l]。i = iabs(i1)。if(i == j)// 支路左節(jié)點(diǎn)號。// 變壓器支路有一節(jié)點(diǎn)號為負(fù)值。//排除左、右節(jié)點(diǎn)號相等的情況。x = zx[l]。// zr[],zx[],zyk[]:支路三參數(shù)。gij = r / zf。bij =1/x。gii[i] = gii[i] + gij。gii[j] = gii[j] + gij。if((i!= iabs(izl[l + 1]))||(j!= iabs(jzl[l + 1])))ll++。對角元實(shí)部為gii，虛部為bii，// 非零非對角元實(shí)部為yb[]，虛部為yb[]，列足碼為iy[]。printf2(gii,bii,n)。double r, x, yk, zf, gij, bij。for(i = 1。bii[j] = bii[j]gij / yk。}if((i!= iabs(izl[l + 1]))||(j!= iabs(jzl[l + 1]))){ll++。}elsekk = 1。對角元實(shí)部為gii，虛部為bii，// 非零非對角元實(shí)部為yb[]，虛部為yb[]，列足碼為iy[]。printf2(gii,bii,n)。for(i = 1。fprintf(fp2, “%%%8d”, yg[i], yb[i], iy[i])。coutUd。//分別為0前面有功，1前面無功，2壓降，3中間有功，4中間無功,5變比,6后面有功，7后面無功int ii。data[ii][0]=。data[ii][2]=Ud。data[ii][4]=Qd。data[ii][6]=。ii++。double b。double d。double f。coutcinx。amp。amp。amp。}if(x==0)break。Node N[100]。N[i].Pn=a。data[ii][0]=。data[ii][2]=。data[ii][4]=N[i].Qn。data[ii][6]