From 209832e59313ea4b9c6d0795f3587079bc3e4b0d Mon Sep 17 00:00:00 2001
From: "dmy@lab" <dmy@lab.com>
Date: Thu, 8 Oct 2015 18:10:39 +0800
Subject: [PATCH] =?UTF-8?q?=E5=8D=95=E7=8B=AC=E5=8A=A0=E4=BA=86=E5=9C=B0?=
 =?UTF-8?q?=E7=BA=BF=E7=9A=84=E8=AE=A1=E7=AE=97=E4=BB=A3=E7=A0=81=EF=BC=8C?=
 =?UTF-8?q?=E7=BB=93=E6=9E=9C=E6=AF=94=E8=BE=83=E6=8E=A5=E8=BF=91=E4=BA=86?=
 =?UTF-8?q?=E3=80=82?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit

Signed-off-by: dmy@lab <dmy@lab.com>
---
 GradualMirror.m | 40 +++++++++++++++++++++++++++++-----------
 1 file changed, 29 insertions(+), 11 deletions(-)

diff --git a/GradualMirror.m b/GradualMirror.m
index e1506b7..0fe3eae 100644
--- a/GradualMirror.m
+++ b/GradualMirror.m
@@ -9,18 +9,25 @@ clear
 semi_lineDistance=450;%分裂间距
 semi_lineCount=6;%分裂数
 ConductorX=[-11000,11000];%导线间距
-ConductorY=[22000,22000,];%导线距地高度
+ConductorY=[22000,22000,30000,30000];%导线距地高度
+GroundX=[-11000,11000];%地线坐标
+GroundY=[30000,30000];
 subconductorR=16.8;%子导线半径
 phaseN=2;%相数，单回三相
+groundN=2;%地线数量
 %%
 eslong=8.854187817*10^-12*1000;
+eslong=1;
 %设置电压
 % Volt_=[1100/sqrt(3);1100/sqrt(3)*exp(1j*4/3*pi);1100/sqrt(3)*exp(1j*2/3*pi);];
-Volt_=[800;-800;];
+Volt_=[800;-800;0;0];
 Volt=[];
 for vLoop=1:phaseN
     Volt=[Volt;Volt_(vLoop)*ones(semi_lineCount,1);];
 end
+for vLoop=1:groundN
+    Volt=[Volt;Volt_(vLoop+phaseN);];
+end
 %按分裂数和分裂导线间距布置单相线路导线
 %用极坐标
 arc=2*pi/semi_lineCount;
@@ -28,21 +35,21 @@ arc=2*pi/semi_lineCount;
 R=semi_lineDistance/2/sin(arc/2);
 %计算导线互电位和自电位系数
 subconductorPos=[];
-for phaseLoop=1:phaseN
+for phaseLoop=1:phaseN%导线
     for sC=1:semi_lineCount
         subconductorPos=[subconductorPos;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 grondLoop=1:groundN%地线
+        subconductorPos=[subconductorPos;GroundX(grondLoop)+1j*GroundY(grondLoop)];%移动到子导线中心
+end
 mirrorSubconductorPos=conj(subconductorPos);%获得子导线镜像
 H=diag(imag(subconductorPos));
-r=eye(length(H))*subconductorR;
+r=eye(length(H))*subconductorR;%暂时认为地线和导线半径一样
 matSubconductor=repmat(subconductorPos,1,length(subconductorPos));
 conductor2conductorDistance=abs(matSubconductor-conj(matSubconductor'));
-% conductor2conductorDistance=conductor2conductorDistance-diag(diag(conductor2conductorDistance));
 conductor2MirrorDistance=abs(matSubconductor-repmat(conj(mirrorSubconductorPos'),length(subconductorPos),1));
-% conductor2MirrorDistance=conductor2MirrorDistance-diag(diag(conductor2MirrorDistance));
 P1=1/2/pi/eslong*log(2*H./r);
 P1(isnan(P1))=0;
 P2=1/2/pi/eslong*log(conductor2MirrorDistance./conductor2conductorDistance);
@@ -67,8 +74,15 @@ for phaseLoop=1:phaseN
         end
     end
 end
+
+for groundLoop=1:groundN%地线只用一根导线等效
+    innerMirrorPos=[innerMirrorPos;subconductorPos(groundLoop+semi_lineCount*phaseN)];%先预留一个位置
+    innerMirrorQ=[innerMirrorQ;sum(QRI(groundLoop+semi_lineCount*phaseN))];
+end
+%以下是验证部分
 %选检验导线上一个角度
-vrfRelA=linspace(0,2*pi,200)';%vrf=verify
+verifyPointN=200;
+vrfRelA=linspace(0,2*pi,verifyPointN)';%vrf=verify
 %计算检验点相对于子导线的位置
 vrfRelPos=exp(1j*vrfRelA)*subconductorR;
 %移动坐标，使验证的子导线中心和实际子导线中心重合。
@@ -78,6 +92,9 @@ for phaseLoop=1:phaseN
         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(innerMirrorPos));
 vrf2ConductorDistance=abs(matVrfPos-repmat(conj(innerMirrorPos'),length(vrfPos),1));
@@ -87,12 +104,13 @@ Pij=1/2/pi/eslong*log(vrf2MirrorDistance./vrf2ConductorDistance);
 V=Pij*innerMirrorQ;
 Vvalidation=[];
 for phaseLoop=1:phaseN
-    Vvalidation=[Vvalidation;Volt_(phaseLoop)*ones(semi_lineCount*200,1);];
+    Vvalidation=[Vvalidation;Volt_(phaseLoop)*ones(semi_lineCount*verifyPointN,1);];
+end
+for groundLoop=1:groundN
+    Vvalidation=[Vvalidation;Volt_(groundLoop+phaseN)*ones(verifyPointN,1);];
 end
 error=abs((V-Vvalidation)./Vvalidation);
-%     error(isinf(error))=0;
 error=sum(error)/length(Vvalidation)
-%以下是验证部分
 display('Finished.');
 %计算场强
 ABCy=imag(repmat(innerMirrorPos,1,length(vrfPos)));