114 lines
4.5 KiB
Matlab
114 lines
4.5 KiB
Matlab
clc
|
|
clear
|
|
%% 自适应模拟电荷法
|
|
% [1]. 任巍巍, 孙.A.宗.A., 一种较准确的分裂导线表面场强计算方法. 电网技术, 2006(04): 第92-96页.
|
|
% [2]. 陈习文, 特高压直流输电线路电磁环境的研究, 2012, 北京交通大学.
|
|
%%
|
|
%设置几个参数
|
|
semi_lineDistance=457;%分裂间距
|
|
semi_lineCount=4;%分裂数
|
|
ConductorX=[-13720,0,13720];%导线距地高度
|
|
ConductorY=[20830,20830,20830];%导线间距
|
|
CSM_N=80;%每一个子导线的模拟电荷数
|
|
subconductorR=29.95;%子导线半径
|
|
phaseN=3;%相数，单回三相
|
|
%%
|
|
%设置电压
|
|
Volt_=[765/sqrt(3);765/sqrt(3)*exp(1j*4/3*pi);765/sqrt(3)*exp(1j*2/3*pi);];
|
|
Volt=[];
|
|
for vLoop=1:length(Volt_)
|
|
Volt=[Volt;Volt_(vLoop)*ones(CSM_N*semi_lineCount,1);];
|
|
end
|
|
%按分裂数和分裂导线间距布置单相线路导线
|
|
%用极坐标
|
|
arc=2*pi/semi_lineCount;
|
|
CSM_arc=2*pi/CSM_N;
|
|
%子导线中心到导线中心的距离
|
|
R=semi_lineDistance/2/sin(arc/2);
|
|
%计算模拟电荷的位置
|
|
r1=20;
|
|
error=10000;
|
|
step=1/10;
|
|
maxLoop=round((subconductorR-r1)/step);
|
|
for Loop=1:maxLoop;
|
|
simulationChargePos=ones(CSM_N,1);
|
|
simulationChargeABCPos=[];
|
|
for I=1:CSM_N
|
|
simulationChargePos(I)=exp(1j*((I-1)*CSM_arc+CSM_arc/2))*r1;%逆时针转一个角度
|
|
end
|
|
for phaseLoop=1:phaseN
|
|
for sC=1:semi_lineCount
|
|
simulationChargeABCPos=[simulationChargeABCPos;simulationChargePos+ConductorX(phaseLoop)+1j*ConductorY(phaseLoop)+exp(1j*((sC-1)*arc+arc/2))*R];%移动到子导线中心
|
|
end
|
|
end
|
|
% simulationChargeAPos=simulationChargePos+ConductorX(1)+1j*ConductorY(1);
|
|
% simulationChargeBPos=simulationChargePos+ConductorX(2)+1j*ConductorY(2);
|
|
% simulationChargePos=simulationChargeABCPos;
|
|
%计算电位系数
|
|
H=diag(imag(simulationChargeABCPos));
|
|
r=subconductorR*eye(length(imag(simulationChargeABCPos)));%导线自几何均距
|
|
%导线与导线的距离
|
|
matSimulationChargePos=repmat(simulationChargeABCPos,1,length(simulationChargeABCPos));
|
|
conductor2conductorDistance=matSimulationChargePos-conj(matSimulationChargePos');
|
|
conductor2conductorDistance=abs(conductor2conductorDistance-diag(diag(conductor2conductorDistance)));
|
|
matMirrorChargePos=conj(matSimulationChargePos);%虚部取负号
|
|
conductor2MirrorDistance=matSimulationChargePos-conj(matMirrorChargePos');
|
|
conductor2MirrorDistance=abs(conductor2MirrorDistance-diag(diag(conductor2MirrorDistance)));
|
|
eslong=8.854187817*10;
|
|
P1=1/2/pi/eslong*log(2*H./r);
|
|
P1(isnan(P1))=0;
|
|
P2=1/2/pi/eslong*log(conductor2MirrorDistance./conductor2conductorDistance);
|
|
P2(isnan(P2))=0;
|
|
P=P1+P2;
|
|
%求电荷
|
|
QRI=P\Volt;
|
|
%以下是验证部分
|
|
if error<0.0001
|
|
break;
|
|
end
|
|
%选检验导线上一个角度
|
|
vrfRelA=linspace(0,2*pi,200)';%vrf=verify
|
|
%计算检验点相对于子导线的位置
|
|
vrfRelPos=exp(1j*vrfRelA)*subconductorR;
|
|
%移动坐标，使验证的子导线中心和实际子导线中心重合。
|
|
vrfPos=[];
|
|
for phaseLoop=1:phaseN
|
|
for sC=1:semi_lineCount
|
|
vrfPos=[vrfPos;exp(1j*((sC-1)*arc+arc/2))*R+ConductorX(phaseLoop)+1j*ConductorY(phaseLoop)+vrfRelPos];
|
|
end
|
|
end
|
|
%计算这一点的电位系数
|
|
matVrfPos=repmat(vrfPos,1,length(simulationChargeABCPos));
|
|
vrf2ConductorDistance=abs(matVrfPos-repmat(conj(simulationChargeABCPos'),length(vrfPos),1));
|
|
vrf2MirrorDistance=abs(matVrfPos-repmat(conj(conj(simulationChargeABCPos')),length(vrfPos),1));
|
|
Pij=1/2/pi/eslong*log(vrf2MirrorDistance./vrf2ConductorDistance);
|
|
%计算电压
|
|
V=Pij*QRI;
|
|
Vvalidation=[];
|
|
for phaseLoop=1:phaseN
|
|
Vvalidation=[Vvalidation;Volt_(phaseLoop)*ones(semi_lineCount*200,1);];
|
|
end
|
|
error=sum(abs((V-Vvalidation)./Vvalidation))/length(Vvalidation);
|
|
r1=r1+step;
|
|
end
|
|
display('Finished.');
|
|
if Loop<maxLoop
|
|
display('Converged.');
|
|
end
|
|
display(Loop);
|
|
%计算场强
|
|
ABCy=imag(repmat(simulationChargeABCPos,1,length(vrfPos)));
|
|
ABCx=real(repmat(simulationChargeABCPos,1,length(vrfPos)));
|
|
y=imag(conj(matVrfPos'));
|
|
x=real(conj(matVrfPos'));
|
|
ERy=sum( ( (ABCy-y)./( (ABCy-y).^2+(ABCx-x).^2 )+(ABCy+y)./( (ABCy+y).^2+(ABCx-x).^2 ) ).*repmat(real(QRI),1,length(vrfPos))./2/pi/eslong,1 );
|
|
EIy=sum( ( (ABCy-y)./( (ABCy-y).^2+(ABCx-x).^2 )+(ABCy+y)./( (ABCy+y).^2+(ABCx-x).^2 ) ).*repmat(imag(QRI),1,length(vrfPos))./2/pi/eslong,1 );
|
|
ERx=sum( ( (ABCx-x)./( (ABCy-y).^2+(ABCx-x).^2 )+(ABCx-x)./( (ABCy+y).^2+(ABCx-x).^2 ) ).*repmat(real(QRI),1,length(vrfPos))./2/pi/eslong,1 );
|
|
EIx=sum( ( (ABCx-x)./( (ABCy-y).^2+(ABCx-x).^2 )+(ABCx-x)./( (ABCy+y).^2+(ABCx-x).^2 ) ).*repmat(imag(QRI),1,length(vrfPos))./2/pi/eslong,1 );
|
|
E=sqrt(ERy.^2+EIy.^2+ERx.^2+EIx.^2);
|
|
max(E);
|
|
|
|
scatter(real(simulationChargeABCPos(1:length(simulationChargeABCPos)/1)),imag(simulationChargeABCPos(1:length(simulationChargeABCPos)/1)),[],'r');
|
|
axis equal
|
|
hold on;
|
|
scatter(real(vrfPos),imag(vrfPos),[],'k'); |