dandaotichuanshuxian.m

一个用matlab编写的天线问题

clear all;close all;
f=0;
for a=0.01:0.005:0.05%半径的变化范围
    f=f+1;
    A(f)=a;
    nn=0;
    H=1;
    N=40;
    for B=0:pi/20:(2*pi-pi/20)
        nn=nn+1;
        X(nn)=a*cos(pi*nn/20);
        Y(nn)=H+a*sin(pi*nn/20);
        l=sqrt(2*a^2-2*a^2*cos(pi/20));
        dx(nn)=l*cos(abs(pi-B-(pi-pi/20)/2));
        dy(nn)=l*sin(abs(pi-B-(pi-pi/20)/2));
    end
    for nn=1:1:N-1
        x(nn)=0.5*(X(nn)+X(nn+1));
        y(nn)=0.5*(Y(nn)+Y(nn+1));
    end
    x(N)=0.5*(X(N-1)+X(1));
    y(N)=0.5*(Y(N-1)+Y(1));
    t1=-0.5,t2=0.5;
    for mm=1:1:40
        for nn=1:1:40
            a1=(x(mm)-x(nn))^2+(y(mm)+y(nn))^2;
            b1=-2*(x(mm)-x(nn))*dx(nn)+2*(y(mm)+y(nn))*dy(nn);
            c1=dx(nn)^2+dy(nn)^2;
            a2=(x(mm)-x(nn))^2+(y(mm)-y(nn))^2;
            b2=-2*(x(mm)-x(nn))*dx(nn)-2*(y(mm)-y(nn))*dy(nn);
            c2=c1;
            d1=b1^2-4*a1*c1;
            d2=b2^2-4*a2*c2;
            if d1>0
                s11=(t1+b1/2*c1)*log(a1+b1*t1+c1*t1^2)-2*t1+(sqrt(d1)/c1)*atanh((2*c1*t1+b1)/sqrt(d1));
                s12=(t2+b1/2*c1)*log(a1+b1*t2+c1*t2^2)-2*t2+(sqrt(d1)/c1)*atanh((2*c1*t2+b1)/sqrt(d1));
            else
                s11=(t1+b1/2*c1)*log(a1+b1*t1+c1*t1^2)-2*t1+(sqrt(-d1)/c1)*atan((2*c1*t1+b1)/sqrt(-d1));
                s12=(t2+b1/2*c1)*log(a1+b1*t2+c1*t2^2)-2*t2+(sqrt(-d1)/c1)*atan((2*c1*t2+b1)/sqrt(-d1));
            end
            s1=s12-s11;
            if d2>0
                s21=(t1+b2/2*c2)*log(a2+b2*t1+c2*t1^2)-2*t1+(sqrt(d2)/c2)*atanh((2*c2*t1+b2)/sqrt(d2));
                s22=(t2+b2/2*c2)*log(a2+b2*t2+c2*t2^2)-2*t2+(sqrt(d2)/c2)*atanh((2*c2*t2+b2)/sqrt(d2));
            else if d2==0
                    s21=(t1+b2/2*c2)*log(a2+b2*t1+c2*t1^2)-2*t1;
                    s22=(t2+b2/2*c2)*log(a2+b2*t2+c2*t2^2)-2*t2;
                else
                    s21=(t1+b2/2*c2)*log(a2+b2*t1+c2*t1^2)-2*t1+(sqrt(-d2)/c2)*atan((2*c2*t1+b2)/sqrt(-d2));
                    s22=(t2+b2/2*c2)*log(a2+b2*t2+c2*t2^2)-2*t2+(sqrt(-d2)/c2)*atan((2*c2*t2+b2)/sqrt(-d2));
                end
            end
            s2=s22-s21;
            s(mm,nn)=l*s1-l*s2;
        end
    end
    for nn=1:1:40
        dL(nn)=l;
        b(nn)=1;
    end
    E=8.854e-12;
    B0=4*pi*E*b;
    C1(f)=dL*(inv(s))*B0';
    Y1=pi*4e-7;
    L1(f)=Y1*E/C1(f);
    Zc1(f)=sqrt(L1(f)/C1(f));
    d=2*a;
    D=2*H;
    C2(f)=2*pi*E/log((D+sqrt(D^2-d^2))/d);
    Y2=pi*4e-7;
    L2(f)=Y2*E/C2(f);
    Zc2(f)=sqrt(L2(f)/C2(f));
end
plot(A,C1,A,C2,'r');%画电容的关系图
xlabel('a(mm)');
ylabel('C(F/mm)');
legend('矩量法结果','解析法结果');
grid on;
figure;
hold on;
plot(A,L1,A,L2,'r');%画电容的关系图
xlabel('a(mm)');
ylabel('L(H/mm)');
legend('矩量法结果','解析法结果');
grid on;
figure;
hold on;
plot(A,Zc1,A,Zc2,'r');%画特性阻抗关系图
xlabel('a(mm)');
ylabel('Zc(Ω)');
grid on;
legend('矩量法结果','解析法结果');