clc clear %% 自适应模拟电荷法 % [1]. 任巍巍, 孙.A.宗.A., 一种较准确的分裂导线表面场强计算方法. 电网技术, 2006(04): 第92-96页. % [2]. 陈习文, 特高压直流输电线路电磁环境的研究, 2012, 北京交通大学. % 采用文献的结果： 蒋兴良, 胡.A.舒.A., Analysis of Conductors＇ Surface Electric Field of UHVDC Transmission Lines Based on Optimized Charge Simulation Method. 高电压技术, 2008(12): p. 2547-2551. %% %设置几个参数 semi_lineDistance=450;%分裂间距 semi_lineCount=6;%分裂数 ConductorX=[-11000,11000];%导线间距 ConductorY=[22000,22000,];%导线距地高度 GroundX=[-11000,11000];%地线坐标 GroundY=[30000,30000]; CSM_N=50;%每一个子导线的模拟电荷数 subconductorR=16.8;%子导线半径 phaseN=2;%相数，单回三相 groundN=2;%地线数量 %% %设置电压 % Volt_=[1100/sqrt(3);1100/sqrt(3)*exp(1j*4/3*pi);1100/sqrt(3)*exp(1j*2/3*pi);]; Volt_=[800;-800;0;0]; Volt=[]; for vLoop=1:phaseN Volt=[Volt;Volt_(vLoop)*ones(CSM_N*semi_lineCount,1);]; end for vLoop=1:groundN Volt=[Volt;Volt_(vLoop+phaseN)*ones(CSM_N,1);]; end %按分裂数和分裂导线间距布置单相线路导线 %用极坐标 arc=2*pi/semi_lineCount; CSM_arc=2*pi/CSM_N; %子导线中心到导线中心的距离 R=semi_lineDistance/2/sin(arc/2); %计算模拟电荷的位置 r1=8; error=10000; step=1/10; maxLoop=round((subconductorR-r1)/step); for Loop=1:maxLoop; simulationChargePos=ones(CSM_N,1); simulationChargeABCPos=[]; matchPos=[]; 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];%移动到子导线中心 %同时计算匹配点的位置 matchPos=[matchPos;simulationChargePos/r1*subconductorR+ConductorX(phaseLoop)+1j*ConductorY(phaseLoop)+exp(1j*((sC-1)*arc+arc/2))*R]; end end %地线的 for groundLoop=1:groundN %simulationChargeABCPos=[simulationChargeABCPos;simulationChargePos+GroundX(phaseLoop)+1j*GroundY(phaseLoop)];%移动到子导线中心 simulationChargeABCPos=[simulationChargeABCPos;simulationChargePos+GroundX(groundLoop)+1j*GroundY(groundLoop)];%移动到子导线中心 %同时计算匹配点的位置 matchPos=[matchPos;simulationChargePos/r1*subconductorR+GroundX(groundLoop)+1j*GroundY(groundLoop)]; % matchPos=[matchPos;GroundX(groundLoop)+1j*GroundY(groundLoop)]; 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)); matMatchPos=repmat(conj(matchPos'),length(simulationChargeABCPos),1); CMS2MatchPointDistance=abs(matSimulationChargePos-matMatchPos); % conductor2conductorDistance=abs(conductor2conductorDistance-diag(diag(conductor2conductorDistance))); matMirrorChargePos=conj(matSimulationChargePos);%虚部取负号 mirrorCharge2MatchPointDistance=abs(matMirrorChargePos-matMatchPos); % conductor2MirrorDistance=abs(conductor2MirrorDistance-diag(diag(conductor2MirrorDistance))); eslong=8.854187817*10^-12*1000; eslong=1; % P1=1/2/pi/eslong*log(2*H./r); % P1(isnan(P1))=0; P2=1/2/pi/eslong*log(mirrorCharge2MatchPointDistance./CMS2MatchPointDistance); % P2(isnan(P2))=0; % P=P1+P2; P=P2; %求电荷 QRI=P\Volt; %选检验导线上一个角度 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 for groundLoop=1:groundN vrfPos=[vrfPos;GroundX(groundLoop)+1j*GroundY(groundLoop)+vrfRelPos]; 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 for groundLoop=1:groundN Vvalidation=[Vvalidation;Volt_(groundLoop+phaseN)*ones(200,1);]; end error=abs((V-Vvalidation)./Vvalidation); error=[error(1:phaseN*semi_lineCount*200);abs(V(phaseN*semi_lineCount*200+1:end)-Vvalidation(phaseN*semi_lineCount*200+1:end))]; % error(isinf(error))=0; error=sum(error)/length(Vvalidation) %以下是验证部分 if error<0.01 break; end r1=r1+1*step; end display('Finished.'); if Loop