改进节点法
#include
#include
#include
//定义复数的实虚部
typedef struct
{
float fl;
float re;
}cplx;
//定义局部变量
int n;
float w;
cplx cmplx(a,b);
float a,b;
{
cplx p;
p.fl=a;
p.re=b;
return p;
}
//复数加法运算的函数
cplx add(a,b)
cplx a,b;
{
cplx p;
p.re=a.re+b.re;
p.fl=a.fl+b.fl;
return p;
}
//复数减法的运算函数
cplx sub(a,b)
cplx a,b;
{
cplx p;
p.re=a.re-b.re;
p.fl=a.fl-b.fl;
return p;
}
//复数乘法运算函数
cplx time(a,b);
cplx a,b;
{
cplx p;
p.re=a.fl*b.re+a.re*b.fl;
p.fl=a.fl*b.fl-a.re*b.re;
return p;
}
//复数除法运算函数
cplx divv(a,b)
cplx a,b;
{
cplx p;
float q;
q=b.re*b.re+b.fl*b.fl;
if(q)
{
q=1/q;
p.re=(a.re*b.fl-a.fl*b.re)*q;
p.fl=(a.fl*b.fl-a.re*b.re)*q;
return p;
}
}
//复数取模运算的函数
float sabs(a)
cplx a;
{
float p;
p=sqrt(a.fl*a.fl+a.re*a.re);
return p;
}
//复数求幅角运算的函数
float ang(a)
cplx a;
{
float p;
if(a.fl==0)
p=90.00;
else
p=atan2(a.re,a.fl);
p=p*180/3.1415926;
return p;
}
//主函数
void main()
{
int input();
float gauss();
void conclus();
int k,i,j,ny;
float deta;
cplx gna[20][20],mv[20].vn[20],amps[20];
cplx b={0,0};
//清零
for(i=0;i
{
jn[i]=b;
mv[i]=b;
vn[i]=b;
for(j=0;j
{
gn[i][j]=b;
gna[i][j]=b;
}
}
//数据输入，形成【ＧＮ］［ＪＮ］［ＧＮＡ］
ny=input(mv);
for(i=1;i
{
for(j=1;j
gna[i][j]=gn[i][j];
gna[i][ny+1]=jn[i];
}
//求解节点方程，结果输出
deta=gauss(gna,vn,ny);
if(deta
printf("\n NO ANSWER!");
else
conclus(vn,mv,ny);
}
//结果输出函数
void conclus(vn,mv,ny)
cplx vn[20];
int mv[20],ny;
{
int i,j;
printf("\n*********节点电压*********\n");
for(i=1;i
printf("VN[%d]=%7.2f+j%7.2f\n",i,vn[i].fl,vn[i].re);
printf("\n*******电压源支路电流*******\n");
for(i=n+1;i
printf("IS[%d]=%7.2f+j%7.2f\n",mv[i],vn[i].fl,vn[i].re);
printf("\n***********时域形式***********\n");
for(i=1;i
printf("vn[%d]=%7.2fcos(%5.0ft+%7.2f)\n",i,sabs(vn[i]),w,ang(vn[i])); for(i=n+1;i
printf("IS[%d]=%7.2fcos(%5.0ft+%7.2f)\n",mv[i],sabs(vn[i]),w,ang(vn[i])); }
//数据输入，形成矩阵[GN][JN]的函数
int input(mv)
int mv[20];
{
void gj();//处理电导，电流源支路函数
void es();//处理电压源支路函数
void vcvs();//处理电压控制电压源支路函数
void vccs();//处理电压控制电流源支路函数
void ccvs();//处理电流控制电压源支路函数
void cccs();//处理电流控制电压源支路函数
void hllm();//处理互感支路函数
int nt,nb,i,j,ns[20],ic[20],ny,nf[20],me[20],a[20][20],jr[20],lcr[20]; float dx[20],dxc[20];
//数组清零
for(i=1;i
{
jr[i]=0;
ns[i]=0;
nf[i]=0;
dx[i]=0;
dxc[i]=0;
me[i]=0;
lcr[i]=0;
mv[i]=0;
ic[i]=0;
for(j=1;j
a[i][j]=0;
}
//输入节点数NT ，支路数NB ，角频率W
scanf("%d",&nt);
scanf("%d",&nb);
scanf("%f",&w);
n=nt-1;
ny=n;
//按 1. 电导支路 2. 电流源支路 3. 电压源支路 4.vcvs 支路 // 5.vccs 支路 6. ccvs支路 7.cccs 支路 8. 互感支路 //输出支路类型
for(i=1;i
{
scanf("%d",&jr[i]);
//调用处理各种类型支路函数，形成矩阵[GN][JN]
swatch(jr[i])
{
case 1:
case 2:gj(i,jr,lcr,ns,nf,dx,dxc,w);
break;
case 3:ny++;
es(i,ns,nf,dx,dxc,ny,me);
mv[ny]=i;
break;
case 4:ny++;
vcvs(i,ns,nf,dx,ny,ic,me);
mv[ny]=i;
break;
case 5:vccs(i,ns,nf,dx,ic);
break;
case 6:ny++;
ccvs(i,ns,nf,dx,ny,ic,me);
mv[ny]=i;
break;
case 7:cccs(i,ns,nf,dx,ic,me);
break;
case 8: hllm(i,ns,nf,dx,dxc,ic,w);
break;
defaut:printf("error");
}
}
for(i=1,i
{
a[ns[i]][i]=1;
a[nf[i]][i]=-1;
}
//打印各支路参数及[A][GN|JN]
printf("\nNT=%d NB=%d NY=%d W=%5.0f\n",nt,nb,ny,w); printf(" I JR NS NF DX ICR DXC IC\n"); for(i=1;i
printf("%4d%5d%5d%5d%8.2f",i,jr[i],ns[i],nf[i],dx[i]); printf("%4d%8.2f%5d\n",lcr[i],dxc[i],ic[i]);
printf("关联矩阵\n");
for(i=1;i
{
for(j=1;j
printf("%4d",a[i][j]);
printf("\n");
}
printf("[GN|JN]%d\n",ny);
for(i=1;i
{
for(j=1;j
printf("%4.1f+j%4.1f",gn[i][j].fl,gn[i][j].re); printf("%4.1f+j%4.1f",jn[i].fl,jn[i].re);
}
return ny;
}
//处理阻抗，电流源支路函数
void gj(i,jr,lcr,ns,nf,dx,dxc,w)
float dxc[20],dx[20],w;
int i,jr[20],nf[20],lcr[20];
{
int j,k;
float x;
cplx a,b,t;
a.fl=1;
a.re=0;
if(jr[i]==1) //阻抗支路
{
//输入起点，止点，参数及元件类型（R ：１，L ：２，C ：３） scanf("%d%d%f%d",&ns[i],&nf[i],&dx[i],&lcr[i]); j=ns[i];
k=nf[i];
if(lcr[i]==1) //电阻
{
b.fl=dx[i];
b.re=0;
}
if(lcr==2) //电感
{
b.fl=0;
b.re=dx[i]*w;
}
if(lcr==3 //电容
{
b.fl=0;
b.re=-1/(dx[i]*w);
}
t=divv(a,b);
if((j!=0)&&(k!=0))
{
gn[j][j]=add(gn[j][j],t); gn[k][k]=add(gn[k][k],t); gn[k][j]=sub(gn[k][j],t); gn[j][k]=sub(gn[j][k],t); }
else if(k==0)
gn[j][j]=add(gn[j][j],t);
else
gn[k][k]=add(gn[k][k],t);
}
else //电流源支路
{
//输入起点，止点与参数及初相角
scanf("%d%d%f%f",ns[i],nf[i],dx[i],dxc[i]); x=dxc[i]*3.14159/180;
b.fl=dx[i]*cos(x);
b.re=dx[i]*sin(x);
j=ns[i];
k=nf[i];
if(j)
jn[j]=sub(jn[j],b);
if(k)
jn[k]=add(jn[k],b);
}
}
//处理电压源支路的函数
void es(i,ns,nf,dx,dxc,ny,me)
int i,ny,nf[20],ns[20],me[20];
float dx[20],dxc[20];
{
int i,j;
cplx b;
cplx gg={1,0};
cplx hh={-1,0};
float x;
//输入电压源的起点，止点与参数及初相角
scanf("%d%d%f%f",&ns[i],&nf[i],&dx[i],&dxc[i]); j=ns[i];
k=nf[i];
me[i]=ny;
if(j!=0)
{
gn[j][ny]=gg;
}
} if(k!=0) { gn[k][ny]=hh; gn[ny][k]=hh; } x=dxc[i]*3.1415926/180; b.fl=dx[i]*cos(x); b.re=dx[i]*sin(x); jn[ny]=b;
//处理VCVS 支路函数
void vcvs(i,ns,nf,dx,ny,ic,me)
int i,ns[20],nf[20],ny,ic[20],me[20];
float dx[20];
{
int j,k,jc,kc;
cplx gg={1,0};
cplx hh={-1,0};
cplx b;
//输入VCVS 起点，止点，参数与控制支路
scanf("%d%d%f%f",&ns[i],&nf[i],&dx[i],&ic[i]); j=ns[i];
k=nf[i];
jc=ns[ic[i]];
kc=nf[ic[i]];
me[i]=ny;
if(j!=0)
{
gn[j][ny]=gg;
gn[ny][j]=gg;
}
if(k!=0)
{
gn[k][ny]=hh;
gn[ny][k]=hh;
}
b.fl=dx[i];
b.re=0;
if(jc!=0)
gn[ny][jc]=sub(gn[ny][jc],b);
if(kc!=0)
}
//处理VCCS 支路函数
void vccs(i,ns,nf,dx,ic)
int i,ns[20],nf[20];
float dx[20];
int ic[20];
{
int j,k,jc,kc;
cplx b;
//输入VCCS 起点止点参数与控制支路 scanf("%d%d%f%d",&ns[i],&nf[i],&dx[i],&ic[i]); j=ns[i]; k=nf[i]; jc=ns[ic[i]]; kc=nf[ic[i]]; b.fl=dx[i]; b.re=0; if(j==0) { if(jc!=0) gn[k][jc]=sub(gn[k][jc],b); if(kc!=0) gn[k][kc]=add(gn[k][kc],b); } else if(k==0) { if(jc!=0) gn[j][jc]=add(gn[j][jc],b); if(kc!=0) gn[j][kc]=sub(gn[j][kc],b); } else { if(jc!=0) { gn[j][jc]=sub(gn[j][kc],b); gn[k][kc]=add(gn[k][kc],b); } if(kc!=0) {
gn[j][kc]=sub(gn[j][kc],b); gn[k][kc]=add(gn[k][kc],b); }
}
}
//处理CCVS 支路函数
void cccs(i,ns,nf,dx,ic,me)
int i,ns[20],nf[20],ic[20],me[20];
float dx[20];
{
int j,k,kc;
cplx gg;
cplx hh;
cplx b;
//输入CCVS 的起点止点，参数控制支路序号 scanf("%d%d%f%d",ns[i],nf[i],dx[i],ic[i]); j=ns[i];
k=nf[i];
kc=me[ic[i]];
me[i]=ny;
gg.fl=1;
gg.re=0;
hh.fl=-1;
hh.re=0;
if(j!=0)
{
gn[j][ny]=gg;
gn[ny][j]=gg;
}
if(k!=0)
{
gn[k][ny]=hh;
gn[ny][k]=hh;
}
b.fl=dx[i];
b.re=0;
gn[ny][kc]=sub(gn[ny][ke],b);
}
//处理CCCS 支路函数
void cccs(i,ns,nf,dx,ic,me)
int i,ns[20],nf[20],ic[20],me[20];
float dx[20];
{
int j,k,kc;
cplx b;
//输入CCCS 的起点，止点，参数及控制支路序号
scanf("%d%d%f%d",ns[i],nf[i],dx[i],ic[i]);
j=ns[i];
k=nf[i];
kc=me[ic[i]];
b.dl=dx[i];
b.re=0;
if(j)
gn[j][kc]=add(gn[j][kc],b);
if(k)
gn[k][kc]=sub(gn[k][kc],b);
}
//处理互感支路函数
void hllm(i,ns,nf,dx,dxc,ic,w)
int i,ns[20],nf[20],ic[20];
float dx[20],dxc[20],w;
{
int j,k,jc,kc;
float dy,d,p,q,pq;
cplx ga,gb,gc,ge,gf;
//输入互感支路的起点，止点，参数，互感系数，控制支路 scanf("%d%d%f%f%d",ns[i],nf[i],dx[i],dxc[i],ic[i]);
if(dx[ic[i]]!=0)
{
j=ns[i];
k=nf[i];
jc=ns[ic[i]];
kc=nf[ic[i]];
dy=dx[ic[i]];
d=dy*dx[i]-dxc[i]*dxc[i];
p=dx[i]/(w*d);
q=dy/(w*d);
pq=dxc[i]/(w*d);
ga=cmplx(0.0,-p);
gb=cmplx(0.0,p);
gc=cmplx(0.0,-q);
gd=cmplx(0.0,q);
ge=cmplx(0.0,pq);
gf=cmplx(0.0,-pq);
if((j!=0)&&(k!=0))
{
gn[j][j]=add(gn[j][j],gc);
gn[k][k]=add(gn[k][k],gc);
gn[j][k]=add(gn[j][k],gd);
gn[k][j]=add(gn[k][j],gd);
}
else if(k==0)
gn[j][j]=add(gn[j][j],gc);
else
gn[k][k]=add(gn[k][k],gc);
if((jc!=0)&&(kc!=0))
{
gn[jc][jc]=add(gn[jc][jc],ga); gn[kc][kc]=add(gn[kc][kc],ga); gn[jc][kc]=add(gn[jc][kc],gb); gn[kc][jc]=add(gn[kc][jc],gb); }
else if(kc==0)
gn[jc][jc]=add(gn[jc][jc],ga);
else
gn[kc][kc]=add(gn[kc][kc],ga);
if((j!=0)&&(jc!=0)
{
gn[j][jc]=add(gn[j][jc],ge); gn[jc][j]=add(gn[jc][j],ge); }
if((k!=0)&&(kc!=0)
{
gn[k][kc]=add(gn[k][kc],ge); gn[kc][k]=add(gn[kc][k],ge); }
if((j!=0)&&(kc!=0)
{
gn[j][kc]=add(gn[j][kc],gf); gn[kc][j]=add(gn[kc][j],gf); }
if((k!=0)&&(jc!=0)
{
gn[k][jc]=add(gn[k][jc],gf); gn[jc][k]=add(gn[jc][k],gf); }
}
}
//求解方程的GAUSS 函数
float gauss(ra,vn,ny)
cplx ra[20][20],vn[20];
int ny;
{
float zu();
int i,j,k,np;
float eps,deta;
cplx alfa[1],beta;
eps=1e-10;
np=ny+1;
for(i=1;i
{
deta=1.0;
deta=zu(alfa,eps,ny,np,i,ra); if(deta
for(j=i;j
ra[i][j]=divv(ra[i][j],alfa[0]); for(k=1;k
if(k!=i) {
} }
}
for(i=1;i
vn[i]=ra[i][np];
return deta;
}
//选列主元函数
float zu(alfa,eps,ny,np,i,ra)
cplx alfa[1],ra[20][20];
float eps;
int i,ny,np;
{
cplx p;
float q;
int ir,m,ic;
q=sabs(ra[i][i]);
m=i;
for(ir=i;ir
{
if(sabs(ra[ir][i])
}
if(q
return 0.0;
alfa[0]=ra[m][i];
if(m==i)
return 1.0;
for(ic=i;ic
{
}
return 1.0;
}
p=ra[m][ic]; ra[m][ic]=ra[i][ic]; ra[i][ic]=p;