電磁數值計算作業(yè)良

1、作業(yè) 13：雙位法分析半波陣子1010020822楊文良已知：已知條件：半波振子天線長為 L ，半徑為a ，且天線長度與波長的關l 求解：電流分布；方向圖；<<L, a <<L = 0.5l , a結構示意圖,圖 1 對稱陣子結構圖1：電流分布圖：2：E 面方向圖：3003：H 面方向圖：4：立體方向性圖：程序%矩量法計算對稱陣子輸入阻抗function Zin = dipole_half(L,a,lambda,N)%L = input('對稱陣子的長度：')%a = input('對稱陣子的半徑：')%lambda = input(&#

2、39;工作波長：')%N U0E0=input('對稱陣子的分段數：') (4*pi)*10-7;(1/(36*pi)*(10-9);%磁導率%介電常數K = 2*pi/lambda;w = K/(U0*E0)0.5);DeltaL=L/(N+1); u=ones(1,N);V=zeros(N,1);V(N+1)/2,1)=1;%xyzlable1 xyzlable2xyzlable3=zeros(N,3);zeros(N,3);zeros(N,3);%1,2,3,4,.%1-,2-,3-.%1+,2+,3+.xyzlable1 xyzlable2xyzlable3xy

3、zlable1(:,1)=0; xyzlable1(:,2)=0;xyzlable1(:,3)=-(N-1)/2*DeltaL:DeltaL:(N-1)/2*DeltaL; xyzlable2(:,1)=0;xyzlable2(:,2)=0;xyzlable2(:,3)=-(N-1)/2*DeltaL-0.5*DeltaL:DeltaL:(N-1)/2*DeltaL-0.5*DeltaL; xyzlable3(:,1)=0;xyzlable3(:,2)=0;xyzlable3(:,3)=-(N-1)/2*DeltaL+0.5*DeltaL:DeltaL:(N-1)/2*DeltaL+0.5*D

4、eltaL;psi1 psi2 psi3 psi4 psi5for=ones(N,N);ones(N,N);ones(N,N);ones(N,N);ones(N,N);%計算 psi1 矩陣%計算 psi2 矩陣%計算 psi3 矩陣%計算 psi4 矩陣%計算 psi5 矩陣k = 1:Nfor m = 1:Nifm = kpsi1(k,m) = 1/(2*pi*DeltaL)*log(DeltaL/a)-j*K/(4*pi); psi2(k,m) = 1/(2*pi*DeltaL)*log(DeltaL/a)-j*K/(4*pi); psi3(k,m) = 1/(2*pi*DeltaL)*

5、log(DeltaL/a)-j*K/(4*pi);elseR1(k,m)=(xyzlable1(k,1)-xyzlable1(m,1)2+(xyzlable1(k,2)-xyzlable1(m,2)2+(x yzlable1(k,3)-xyzlable1(m,3)2)0.5;psi1(k,m) =exp(-j*K*R1(k,m)/(4*pi*R1(k,m);R2(k,m)=(xyzlable2(k,1)-xyzlable2(m,1)2+(xyzlable2(k,2)-xyzlable2(m,2)2+(x yzlable2(k,3)-xyzlable2(m,3)2)0.5;psi2(k,m) =

6、exp(-j*K*R2(k,m)/(4*pi*R2(k,m);R3(k,m)=(xyzlable3(k,1)-xyzlable3(m,1)2+(xyzlable3(k,2)-xyzlable3(m,2)2+(x yzlable3(k,3)-xyzlable3(m,3)2)0.5;psi3(k,m) =exp(-j*K*R3(k,m)/(4*pi*R3(k,m);endendendfor k = 1:Nfor m = 1:Nif(k-m)=-1psi4(k,m) = 1/(2*pi*DeltaL)*log(DeltaL/a)-j*K/(4*pi); elseR4(k,m)=(xyzlable3(

7、k,1)-xyzlable2(m,1)2+(xyzlable3(k,2)-xyzlable2(m,2)2+(x yzlable3(k,3)-xyzlable2(m,3)2)0.5;psi4(k,m) =exp(-j*K*R4(k,m)/(4*pi*R4(k,m);endendendfor k = 1:Nfor m = 1:Nif(k-m)=1psi5(k,m) = 1/(2*pi*DeltaL)*log(DeltaL/a)-j*K/(4*pi);elseR5(k,m)=(xyzlable2(k,1)-xyzlable3(m,1)2+(xyzlable2(k,2)-xyzlable3(m,2)2

8、+(x yzlable2(k,3)-xyzlable3(m,3)2)0.5;psi5(k,m) =exp(-j*K*R5(k,m)/(4*pi*R5(k,m);endendendZ = ones(N); for k = 1:Nfor m = 1:N%計算 Z 矩陣Z(k,m)=j*w*U0*DeltaL2*psi1(k,m)+1/(j*w*E0)*(psi2(k,m)+psi3(k,m)-psi4(k,m)-psi5 (k,m);endend Ii=inv(Z)*V;Zin=1/Ii(N+1)/2);%對稱陣子輸入阻抗%矩量法計算對稱陣子參數 雙位法 E 面方向圖function Zin =

9、dipole_half_E(L,a,lambda,N)%L = input('對稱陣子的長度：')%a = input('對稱陣子的半徑：')%lambda = input('工作波長：')%NU0 E0=input('對稱陣子的分段數：')(4*pi)*10-7; (1/(36*pi)*(10-9);%磁導率%介電常數K = 2*pi/lambda;w = K/(U0*E0)0.5);DeltaL=L/(N+1); u=ones(1,N);V=zeros(N,1);V(N+1)/2,1)=1;xyzlable1 xyzlable

10、2xyzlable3=zeros(N,3);zeros(N,3);zeros(N,3);%1,2,3,4,.%1-,2-,3-.%1+,2+,3+.xyzlable1 xyzlable2xyzlable3xyzlable1(:,1)=0; xyzlable1(:,2)=0;xyzlable1(:,3)=-(N-1)/2*DeltaL:DeltaL:(N-1)/2*DeltaL; xyzlable2(:,1)=0;xyzlable2(:,2)=0;xyzlable2(:,3)=-(N-1)/2*DeltaL-0.5*DeltaL:DeltaL:(N-1)/2*DeltaL-0.5*DeltaL;

11、 xyzlable3(:,1)=0;xyzlable3(:,2)=0;xyzlable3(:,3)=-(N-1)/2*DeltaL+0.5*DeltaL:DeltaL:(N-1)/2*DeltaL+0.5*DeltaL;psi1 psi2 psi3 psi4 psi5for=ones(N,N);ones(N,N);ones(N,N);ones(N,N);ones(N,N);%計算 psi1 矩陣%計算 psi2 矩陣%計算 psi3 矩陣%計算 psi4 矩陣%計算 psi5 矩陣k = 1:Nfor m = 1:Nifm = kpsi1(k,m) = 1/(2*pi*DeltaL)*log(

12、DeltaL/a)-j*K/(4*pi); psi2(k,m) = 1/(2*pi*DeltaL)*log(DeltaL/a)-j*K/(4*pi); psi3(k,m) = 1/(2*pi*DeltaL)*log(DeltaL/a)-j*K/(4*pi);elseR1(k,m)=(xyzlable1(k,1)-xyzlable1(m,1)2+(xyzlable1(k,2)-xyzlable1(m,2)2+(x yzlable1(k,3)-xyzlable1(m,3)2)0.5;psi1(k,m) =exp(-j*K*R1(k,m)/(4*pi*R1(k,m);R2(k,m)=(xyzlabl

13、e2(k,1)-xyzlable2(m,1)2+(xyzlable2(k,2)-xyzlable2(m,2)2+(x yzlable2(k,3)-xyzlable2(m,3)2)0.5;psi2(k,m) =exp(-j*K*R2(k,m)/(4*pi*R2(k,m);R3(k,m)=(xyzlable3(k,1)-xyzlable3(m,1)2+(xyzlable3(k,2)-xyzlable3(m,2)2+(x yzlable3(k,3)-xyzlable3(m,3)2)0.5;psi3(k,m) =exp(-j*K*R3(k,m)/(4*pi*R3(k,m);endendendfor k

14、 = 1:Nfor m = 1:Nif(k-m)=-1psi4(k,m) = 1/(2*pi*DeltaL)*log(DeltaL/a)-j*K/(4*pi); elseR4(k,m)=(xyzlable3(k,1)-xyzlable2(m,1)2+(xyzlable3(k,2)-xyzlable2(m,2)2+(x yzlable3(k,3)-xyzlable2(m,3)2)0.5;psi4(k,m) =exp(-j*K*R4(k,m)/(4*pi*R4(k,m);endendendfor k = 1:Nfor m = 1:Nif(k-m)=1psi5(k,m) = 1/(2*pi*Delt

15、aL)*log(DeltaL/a)-j*K/(4*pi);elseR5(k,m)=(xyzlable2(k,1)-xyzlable3(m,1)2+(xyzlable2(k,2)-xyzlable3(m,2)2+(x yzlable2(k,3)-xyzlable3(m,3)2)0.5;psi5(k,m) =exp(-j*K*R5(k,m)/(4*pi*R5(k,m);endendendZ = ones(N); for k = 1:Nfor m = 1:N%計算 Z 矩陣Z(k,m)=j*w*U0*DeltaL2*psi1(k,m)+1/(j*w*E0)*(psi2(k,m)+psi3(k,m)-

16、psi4(k,m)-psi5(k,m);endend Ii=inv(Z)*V;Zin=1/Ii(N+1)/2);%對稱陣子輸入阻抗phi=0;theta=-pi:pi/360:pi;%E 面和 H 面方向圖的描繪index=ones(length(theta),N); for m=1:length(theta)for n=1:Nindex(m,n)=j*K*(xyzlable1(n,1)*sin(theta(m)*cos(phi)+xyzlable1(n,2)*sin(theta( m)*sin(phi)+xyzlable1(n,3)*cos(theta(m);endF(1,m)=DeltaL

17、*(-sin(theta(m)*exp(index(m,:)*Ii;end Fmax=max(abs(F);F1=abs(F)./Fmax x=find(F1<0.709)&(F1>0.705) theta1=-180:0.5:180;y= theta1(x(2)-theta1(x(1) F2=abs(cos(pi/2*cos(theta1*pi/180)./sin(theta1*pi/180);F2(1 length(theta)=0;figure(1)%E 面歸一化方向圖（直角坐標）plot(theta1,F1,'r-',theta1,F2,'b

18、:','linewidth',2)hold onplot(theta1(x(1) theta1(x(1),0plot(theta1(x(2) theta1(x(2),0 cell_string11='fontsize14' cell_string12='it 半功率波瓣寬度text(-90,0.7,cell_string1)F1(x(1),'b-.')F1(x(2),'b-.')= it80circ'legend('雙位法','法',0)title('fontsize

19、14半波陣子 E 面歸一化方向圖') axis(-180,180,0,1.05)xlabel('fontsize14itthetacirc'),ylabel('fontsize14it|E/E_max|')text(0,0.5,'fontsize16phi=0')set(gca,'Xtick',-180,theta1(x(1),-90,theta1(x(2),0, 90,180,'Ytick',0,0.707,1)grid onfigure(2)%極坐標polar(theta,F1,'r-'

20、) hold on polar(theta,F2,'b:')legend('雙位法','法',0)title('fontsize14半波陣子 E 面歸一化方向圖') xlabel('fontsize14itthetacirc') text(0,0.5,'fontsize14it|E/E_max|')%axis(-pi,pi,0,1)%矩量法計算對稱陣子電流分布function Zin = dipole_half_I(L,a,lambda,N)%L = input('對稱陣子的長度：'

21、)%a = input('對稱陣子的半徑：')%lambda = input('工作波長：')%N U0E0=input('對稱陣子的分段數：') (4*pi)*10-7;(1/(36*pi)*(10-9);%磁導率%介電常數K = 2*pi/lambda;w = K/(U0*E0)0.5);DeltaL=L/(N+1); u=ones(1,N);V=zeros(N,1);V(N+1)/2,1)=1;xyzlable1 xyzlable2 xyzlable3xyzlable4=zeros(N,3);zeros(N,3);zeros(N,3);ze

22、ros(N+2,1);%1,2,3,4,.%1-,2-,3-.%1+,2+,3+.xyzlable1 xyzlable2xyzlable3xyzlable1(:,1)=0; xyzlable1(:,2)=0;xyzlable1(:,3)=-(N-1)/2*DeltaL:DeltaL:(N-1)/2*DeltaL; xyzlable2(:,1)=0;xyzlable2(:,2)=0;xyzlable2(:,3)=-(N-1)/2*DeltaL-0.5*DeltaL:DeltaL:(N-1)/2*DeltaL-0.5*DeltaL; xyzlable3(:,1)=0;xyzlable3(:,2)=

23、0;xyzlable3(:,3)=-(N-1)/2*DeltaL+0.5*DeltaL:DeltaL:(N-1)/2*DeltaL+0.5*DeltaL;psi1 psi2psi3=ones(N,N);ones(N,N);ones(N,N);%計算 psi1 矩陣%計算 psi2 矩陣%計算 psi3 矩陣psi4 = ones(N,N); %計算 psi4 矩陣psi5 = ones(N,N); %計算 psi5 矩陣for k = 1:Nfor m = 1:Nif m = kpsi1(k,m) = 1/(2*pi*DeltaL)*log(DeltaL/a)-j*K/(4*pi); psi2

24、(k,m) = 1/(2*pi*DeltaL)*log(DeltaL/a)-j*K/(4*pi); psi3(k,m) = 1/(2*pi*DeltaL)*log(DeltaL/a)-j*K/(4*pi);elseR1(k,m)=(xyzlable1(k,1)-xyzlable1(m,1)2+(xyzlable1(k,2)-xyzlable1(m,2)2+(x yzlable1(k,3)-xyzlable1(m,3)2)0.5;psi1(k,m) =exp(-j*K*R1(k,m)/(4*pi*R1(k,m);R2(k,m)=(xyzlable2(k,1)-xyzlable2(m,1)2+(x

25、yzlable2(k,2)-xyzlable2(m,2)2+(x yzlable2(k,3)-xyzlable2(m,3)2)0.5;psi2(k,m) =exp(-j*K*R2(k,m)/(4*pi*R2(k,m);R3(k,m)=(xyzlable3(k,1)-xyzlable3(m,1)2+(xyzlable3(k,2)-xyzlable3(m,2)2+(x yzlable3(k,3)-xyzlable3(m,3)2)0.5;psi3(k,m) =exp(-j*K*R3(k,m)/(4*pi*R3(k,m);endendendfor k = 1:Nfor m = 1:Nif(k-m)=-

26、1psi4(k,m) = 1/(2*pi*DeltaL)*log(DeltaL/a)-j*K/(4*pi); elseR4(k,m)=(xyzlable3(k,1)-xyzlable2(m,1)2+(xyzlable3(k,2)-xyzlable2(m,2)2+(x yzlable3(k,3)-xyzlable2(m,3)2)0.5;psi4(k,m) =exp(-j*K*R4(k,m)/(4*pi*R4(k,m);endendendfor k = 1:Nfor m = 1:Nif(k-m)=1psi5(k,m) = 1/(2*pi*DeltaL)*log(DeltaL/a)-j*K/(4*p

27、i);elseR5(k,m)=(xyzlable2(k,1)-xyzlable3(m,1)2+(xyzlable2(k,2)-xyzlable3(m,2)2+(xyzlable2(k,3)-xyzlable3(m,3)2)0.5;psi5(k,m) =exp(-j*K*R5(k,m)/(4*pi*R5(k,m);endendendZ = ones(N); for k = 1:Nfor m = 1:N%計算 Z 矩陣Z(k,m)=j*w*U0*DeltaL2*psi1(k,m)+1/(j*w*E0)*(psi2(k,m)+psi3(k,m)-psi4(k,m)-psi5 (k,m);endend

28、 Ii=inv(Z)*V;Zin=1/Ii(N+1)/2);% 電流分布曲線 與正弦電流對比xyzlable4(1 N+2) =-L/2,L/2; xyzlable4(2:N+1)=xyzlable1(:,3); I(1 N+2)=0;I(2:N+1)=abs(Ii)./max(abs(Ii);x = (xyzlable4+L/2)*pi/L;y = sin(x);plot(xyzlable4,I,'b-',xyzlable4,y,'r-.','linewidth',2) xlabel('fontsize14itz'),ylabe

29、l('fontsize14it|I/I_max|') title('fontsize14半波對稱陣子電流分布')axis(-L/2,L/2,0,1)legend('雙位法','正弦電流',0) grid on%輸入阻抗,增益隨半徑頻率特性c=3.0e8f0 = 1.5e8 f=f0*0:0.01:3a = 0.000001*1 10 100 1000 10000 100000;for n=1:length(a) for m=1:length(f)Zin(n,m) = real(dipole_half(1,a(n),c/f(m),3

30、);Xin(n,m) = imag(dipole_half(1,a(n),c/f(m),3);endendform=1:length(f)G(1,m)= dipole_half_G(1,0.01,c/f(m),3,1,pi/2,0)endratio=0:0.01:3;figure(1)plot(ratio,Zin(1,:),'b-',ratio,Zin(2,:),'g-',ratio,Zin(3,:),'r-',ratio,Zin(4,:),'c-',ratio,Zin(5,:),'m-',ratio,Zin(6,:),'k-','linewidth',2)legend('a/L =