矩分析法matlab,用矩量法分析V型天线的代码在matlab2008里老报错-CSDN博客

本帖最后由 msnzest 于 2013-4-21 08:44 编辑

是代码不兼容么？求高手指点个方向，应该怎么才能让它运行起来

clc;

clear all;

close all;

N=40; %每根振子分为N段,N必为偶数 100段的结果比20段的结果接近理论值

lamda=1; %波长

a=0.00001*lamda; %天线的半径

%L=1.5*lamda; %天线长度

L=lamda/2;

%jd=59.25; %两天线之间夹角的一半

jd=90;

deltz=L/N;

e=1e-9/(36*pi); %介电常数epxim

k=2*pi/lamda; %波数

c=2.997924574e8; %光速

omiga=2*pi*c/lamda; %角速度

mui=4*pi*(1e-7); %磁导率

ada=sqrt(mui/e); %波阻抗

pijd=(90-jd)*pi/180;

[zzb,R,Rzf,Rfz,Rff,Rzz,deltz]=zb(N,jd,deltz,L,a);

for m=1:N-1;

for n=1:N-1;

if (m==n) % 1,1 2,2

Psai(m,n)=log(deltz/a)/(2*pi*deltz)-(i*k/(4*pi));

else

Psai(m,n)=exp(-i*k*R(m,n))/(4*pi*R(m,n));

end

if (m==n) % 2-,2-

Psaiff(m,n)=log(deltz/a)/(2*pi*deltz)-(i*k/(4*pi));

else

Psaiff(m,n)=exp(-i*k*Rff(m,n))/(4*pi*Rff(m,n));

end

if (m==n) % 3+,3+

Psaizz(m,n)=log(deltz/a)/(2*pi*deltz)-(i*k/(4*pi));

else

Psaizz(m,n)=exp(-i*k*Rzz(m,n))/(4*pi*Rzz(m,n));

end

if (n==m+1) % 3+,4- 3+=4-

Psaizf(m,n)=log(deltz/a)/(2*pi*deltz)-(i*k/(4*pi));

else

Psaizf(m,n)=exp(-i*k*Rzf(m,n))/(4*pi*Rzf(m,n));

end

if (n==m-1) % 3+,2- 3+=2-

Psaifz(m,n)=log(deltz/a)/(2*pi*deltz)-(i*k/(4*pi));

else

Psaifz(m,n)=exp(-i*k*Rfz(m,n))/(4*pi*Rfz(m,n));

end

for m=1:N-1;

for n=1:N-1;

Z(m,n)=i*omiga*mui*deltz*deltz*Psai(m,n)+(Psaizz(m,n)-Psaizf(m,n)-Psaifz(m,n)+Psaiff(m,n))/(i*omiga*e);

%Z(m,n)=i*omiga*mui*zzb(m,5)*zzb(n,5)*Psai(m,n)+((2*Psai(m,n)-Psaizf(m,n)-Psaifz(m,n))/(i*omiga*e)); %近似认为Psai(m,n)和Psai(m+,n+)和Psai(m-,n-)相等,计算结果基本一样

end

V=zeros(N-1,1); % 电压只有馈电点电压为1，其余点都为0。

V(N/2)=1;

I=Z\V;

Zin=1/I(N/2) % 馈电点的输入阻抗

zb=0:N;

I1=zeros(N+1,1);

for mm=2:N;

I1(mm)=I(mm-1);

end

figure; %各段电流值大小表示图

plot(zb,abs(I1),'r'),title('The current distrubition of The Vee Dipole'),xlabel('the number of segments '),ylabel(' current I (A)');

grid on;

text(30,0.009,'\fontsize{10}L /{\lambda}=0.5');

text(30,0.01,'\fontsize{10}{\gamma}=180^o ');

text(26,0.0115,'Two-potential integral equation');

syms Thita1 Fai;

for nn=1:N-1;

%与每小段相距单位长度(r＝1)上单位电流产生的辐射电场值

%%%%%%%%%%%%%%%% l变换为x+y+z

E_t(1,nn)=i*omiga*mui*exp(-i*k)*(sin(pijd)*cos(Thita1)*sin(Fai)-cos(pijd)*sin(Thita1))*exp(i*k*(zzb(nn,5)*sin(Thita1)*sin(Fai)+zzb(nn,2)*cos(Thita1)))*deltz/(4*pi);

E_f(1,nn)=i*omiga*mui*exp(-i*k)*sin(pijd)*cos(Fai)*exp(i*k*(zzb(nn,5)*sin(Thita1)*sin(Fai)+zzb(nn,2)*cos(Thita1)))*deltz/(4*pi);

E_n(1,nn)=sqrt(E_t(1,nn).*conj(E_t(1,nn)+E_f(1,nn).*conj(E_f(1,nn))));

end; %电场矩阵E(r,Thita,Fai)-- 整个天线在距离单位长度上产生的辐射电场值

E=E_n*I;

for dd=0:0.1:pi;

max=subs(E,{Thita1,Fai},{dd,[0:0.1:2*pi]});

end

% E的最大值,Thita和Fai都取pi/2

D=4*pi*abs(E_max)*abs(E_max)/(120*pi*real(Zin)*abs(I(N/2,1))*abs(I(N/2,1))); %增益系数

G=10*log10(D);

F=abs(E)/E_max; %归一化方向函数

F1=subs(F,{Thita1,Fai},{pi/2,pi/2});