工频电场计算等效电荷法

#include 
#include
#include
#define n 5
#define k 5
double pi=3.1415926;
//gd是导线高度，xd是相间距，dh是导线的垂直距离差 
void dianweixishu（double （*b）[k]double gddouble xddouble dh）
{
double a[n][k];
double h[n]l[n][k]ll[n][k]r[n];
double hls;
int ij;
h[0]=2*（gd+dh）;h[1]=2*gd;h[2]=2*（gd+dh）;
r[0]=0.4;r[1]=0.4;r[2]=0.4;//单位是米
l[0][1]=l[1][0]=l[1][2]=l[2][1]=sqrt（（xd/2）*（xd/2）+dh*dh）;
l[0][2]=l[2][0]=xd;
ll[0][1]=ll[1][0]=ll[1][2]=ll[2][1]=sqrt（（xd/2）*（xd/2）+（2*gd+dh）*（2*gd+dh））;
ll[0][2]=ll[2][0]=sqrt（（2*（gd+dh））*（2*（gd+dh））+xd*xd）;

for（i=0;i<3;i++）
{
	for（j=0;j<3;j++）
{
	
	if （i==j）
	{
		a[i][j]=1.8*pow（1010）*（log（2*h[i]/r[i]））/log（2.718）;
		
		
	}
	else
    {
        a[i][j]=1.8*pow（1010）*（log（ll[i][j]/l[i][j]））/log（2.718）;
		
	}
	
}
        	
		
	
}
    
	hls=a[0][0]*a[1][1]*a[2][2]+a[0][1]*a[1][2]*a[2][0]+a[0][2]*a[1][0]*a[2][1]-a[0][2]*a[1][1]*a[2][0]-a[0][0]*a[2][1]*a[1][2]-a[0][1]*a[1][0]*a[2][2];
    if （hls==0）
    {
		printf（“矩阵不能求逆！数据错误！“）;
    }
	else
	{
		b[0][0]=（a[1][1]*a[2][2]-a[1][2]*a[2][1]）/hls;
		b[0][1]=-（a[0][1]*a[2][2]-a[0][2]*a[2][1]）/hls;
		b[0][2]=（a[0][1]*a[1][2]-a[0][2]*a[1][1]）/hls;
        b[1][0]=-（a[1][0]*a[2][2]-a[1][2]*a[2][0]）/hls;
        b[1][1]=（a[0][0]*a[2][2]-a[0][2]*a[2][0]）/hls;
        b[1][2]=-（a[0][0]*a[1][2]-a[0][2]*a[1][0]）/hls;
        b[2][0]=（a[1][0]*a[2][1]-a[2][0]*a[1][1]）/hls;
		b[2][1]=-（a[0][0]*a[2][1]-a[0][1]*a[2][0]）/hls;
		b[2][2]=（a[0][0]*a[1][1]-a[0][1]*a[1][0]）/hls;
	}
  
  
  	
}
void main（）
{
 double gd;
double xd;
double dh;

double xy; 
double u1[n]u2[n]u3[n]q