%  MOM_for_Antenna.m
%  本函數用于實現用矩量法和海倫積分方程球解對稱天線電流分布
%  編程人    沙威(Wei Sha) 安徽大學(Anhui University) ws108@ahu.edu.cn


%  wave_length_value 波長
%  measurement 物體電尺寸
%  v0 電壓常數
%  count_point 待求點數
%  range 求解的天線歸一化長度
%  a 導體半徑
%  k 波數
%  l 天線長度
%  step 點匹配間距

clear;clc;

range=1/4;
point=3;
wave_length_value=1;
measurement=7.022*(10)^(-3);
v0=1;

k=(2*pi)/wave_length_value;
a=measurement*wave_length_value;
l=range*wave_length_value;
step=l/(point-1);

i_point=1:point;  %  測試點賦值
matrix_wavelength(i_point)=step*(i_point-1);  

z=linspace(-l,l,100);  %  積分離散

%  求解A向量
for i_point=1:point;
  r=((matrix_wavelength(i_point)-z).^2 + a^2).^(1/2);  %  場源距離離散
  g=exp(-j*k*r)./r;  %  格林函數離散
  A(i_point)=trapz(z,cos(k*z).*g);  %  A元素確定
end

%  求解B向量
for i_point=1:point;
  r=((matrix_wavelength(i_point)-z).^2 + a^2).^(1/2);  %  場源距離離散
  g=exp(-j*k*r)./r;  %  格林函數離散
  B(i_point)=trapz(z,sin(2*k*abs(z)).*g);  %  B元素確定
end

%  求解C向量
for i_point=1:point;
   C(i_point)=cos(k*matrix_wavelength(i_point)); %  B元素確定
end

%  阻抗矩陣確定
impedance_matrix=[A.',B.',C.'];

%  電壓矩陣確定
for i_point=1:point;
   voltage_matrix(i_point)=(-j*v0/60)*sin(k*abs(matrix_wavelength(i_point))); %  B元素確定
end

%  求解a1,a2,C_contant;
current=impedance_matrix\voltage_matrix';

%  圖形表示
z_distribute=linspace(0,l,100);
current_function=current(1,1)*sin(k*(l-abs(z_distribute))) +...
                 current(2,1)*sin(2*k*(l-abs(z_distribute)));  %  離散化電流分布

%  電流實虛部
current_re=real(current_function);
current_im=imag(current_function);


%  繪圖
plot(current_re,z_distribute,'r');
hold on;
plot(current_im,z_distribute,'g');

xlabel('current distribution');
ylabel('unitary distance')
title('antenna current distribution plot');
legend('real current','imag current',2);