TE.zip_TE_TE scattering_te波_te波入射_散射

% TE波入射金属椭圆柱体的散射 function TE syms x L S; Z=377; lamda=3; K=2*pi/lamda; a=0.25*lamda; %a,b为椭圆的长轴和短轴 b=1*lamda; N=720; %剖分数目 fai=0; %入射角 Q=(0:2*pi/N:2*pi*(N-1)/N); h=2*pi/N; X=a*cos(Q+h/2); Y=b*sin(Q+h/2); h=2*pi/N; P=zeros(N,N); fun=inline('sqrt((L)^2*(sin(x).^2)+(S)^2*(cos(x).^2))','x','L','S'); %椭圆线长 for n=1:N %通过数值积分计算椭圆每个剖分的线长 x1=(n-1)*h;x2=n*h; %积分上下限 d=quadl(fun,x1,x2,{},{},a,b);%积分 delta(n)=d; end for n=1:N sinn=sign((cos(Q(n)+h/2)))*sqrt(1/((a*tan(Q(n)+h/2)/b)^2+1)); cosn=-1*sign((sin(Q(n)+h/2)))*sqrt(1/(1+(b^2/(a*tan(Q(n)+h/2))^2))); Xn=a*cos(Q(n)+h/2);Yn=b*sin(Q(n)+h/2); for m=1:N Xm=a*cos(Q(m)+h/2);Ym=b*sin(Q(m)+h/2); Rmn=sqrt((Xm-Xn)^2+(Ym-Yn)^2); if m==n P(m,n)=1/2; else P(m,n)=-K*delta(n)*(sinn*(Xm-Xn)/Rmn-cosn*(Ym-Yn)/Rmn)*besselh(1,2,K*Rmn)/(4*i); end end end for n=1:N M(n)=exp(i*K*(X(n)*cos(fai)+Y(n)*sin(fai))); end [L,U]=lu(P); j=U\(L\(M.')); %用直接LU分解求解电流密度 j=j./(M.'); w2=(0:N/2-1); figure(1) plot(w2*2/N,abs(j(1:N/2)),'.'); grid on title('TE波入射金属椭圆体的等效表面电流密度分布图') xlabel('S') ylabel('abs(J/H),current dense/incident magnetic field'); w=(0:pi/180:pi*179/180); out=zeros(1,180); for m=1:180 %计算散射截面（回波宽度） for n=1:N sinn=sign((cos(Q(n)+h/2)))*sqrt(1/((a*tan(Q(n)+h/2)/b)^2+1)); cosn=-1*sign((sin(Q(n)+h/2)))*sqrt(1/(1+(b^2/(a*tan(Q(n)+h/2))^2))); v=exp(i*K*(X(n)*cos(w(m))+Y(n)*sin(w(m)))).*delta(n); g(n)=v*(sinn*cos(w(m))-cosn*sin(w(m))); end G=g*j; T=abs(sum(G)); out(m)=K*T^2/4; end W=(0:1:179); figure(2) plot(W,sqrt(out(180:-1:1)/lamda),'.') %画图，散射截面，0－180度 grid on title('TE波入射金属圆柱体的散射截面分布图') xlabel('Degrees') ylabel('square root of echo width/wavelength')