MOM.rar_MOM.rar_hallen_matlab矩量法_电流分布_矩量法 Hallen

%C***************************************************************** %C MoM %C***************************************************************** %C This is a MATLAB based program using %C %C I. POCKLINGTON'S INTEGRAL EQUATION (8-24) %C II. HALLEN'S INTEGRAL EQUATION (8-27) %C %C That computes: %C %C A. CURRENT DISTRIBUTION %C B. INPUT IMPEDANCE %C C. NORMALIZED AMPLITUDE RADIATION PATTERN %C %C Of a linear symmetrically-excited dipole. %C %C THIS PROGRAM USES PULSE EXPANSION FOR THE ELECTRIC CURRENT MODE %C AND POINT-MATCHING FOR THE ELECTRIC FIELD AT THE CENTER OF EACH %C WIRE SEGMENT. %C %C DELTA-GAP FEED MODEL IS USED IN BOTH FORMULATIONS. IN ADDITION, %C MAGNETIC-FRILL GENERATOR IS AVAILABLE IN THE POCKLINGTON'S %C INTEGRAL EQUATION. %C %C OPTION I. POCKLINGTON'S INTEGRAL EQUATION %C OPTION II. HALLEN'S INTEGRAL EQUATION %C %C %C ** INPUT DATA: %C %C TL = TOTAL DIPOLE LENGTH (IN WAVELENGTHS) %C RA = RADIUS OF THE WIRE (IN WAVELENGTHS) %C NM = TOTAL NUMBER OF SUBSECTIONS (MUST BE AN ODD INTEGER) %C IEX = OPTION TO USE EITHER MAGNETIC-FRILL GENERATOR OR DELTA GAP %C %C IEX = 1 : MAGNETIC-FRILL GENERATOR %C IEX = 2 : DELTA-GAP FEED %C %C ** NOTE: IGNORE INPUT PARAMETER IEX WHEN CHOOSING OPTION II %C (i.e., HALLEN'S FORMULATION) %C ***************************************************************** clear; close; %********************************************************* %*** OUTPUT DEVICE OPTION # : 1 ---> SCREEN % 2 ---> OUTPUT FILE %********************************************************* disp('OUTPUT DEVICE OPTION FOR THE OUTPUT PARAMETERS'); disp(' OPTION (1):SCREEN'); disp(' OPTION (2):OUTPUT FILE'); option_a=0; option_b=0; filename=0; while((option_a~=1)&(option_a~=2)) option_a=input('OUTPUT DEVICE ='); if(option_a==2) filename=input('INPUT THE DESIRED OUTPUT FILENAME <in single quotes> = ','s'); end end %********************************************************* %*** READING OPTION # : 1 ---> POCKLINGTON'S EQUATIONS % 2 ---> HALLEN'S EQUATION %********************************************************* disp('CHOICE OF POCKLINGTON''S OR HALLEN''S EQN.'); disp(' OPTION (1):POCKLINGTON''S EQN.'); disp(' OPTION (2):HALLEN''S EQN.'); option=input('OPTION NUMBER ='); switch option %^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ %*********************************************************************** %*** %POCKLINGTON'S INTEGRAL EQUATIONS *** %*********************************************************************** case 1 %****************************************************** % SOME CONSTANTS AND ARRAYS %****************************************************** nmax=200; nmt=2*nmax-1; cgaw=zeros(nmax); cga=cgaw(1,1:nmax); zmnw=zeros(nmt); zmn=zmnw(1,1:nmt); waw=zeros(nmt); wa=waw(1,1:nmt); etm=zeros(361); etmm=etm(1,1:361); pi=3.14159265; rad=pi/180; beta=2.0*pi; eta=120*pi; %****************************************************** %*** INPUT DATA %****************************************************** nm=input('NUMBER OF SUBDIVISIONS <ODD NUMBER> ='); tl=input('TOTAL DIPOLE LENGTH <WAVELENGTHS> ='); ra=input('RADIUS OF DIPOLE <WAVELENGTHS> ='); %********************************************************* %*** EXITATION OPTION # : 1 ---> MAGNETIC-FRILL % 2 ---> DELTA-GAP %********************************************************** disp('EXCITATION :'); disp(' OPTION (1):MAGNETIC-FRILL'); disp(' OPTION (2):DELTA-GAP'); iex=input('OPTION NUMBER :'); hl=tl/2; nmh=0.5*(nm+1); dz=2*hl/nm; zm=hl-0.5*dz; b=0.5*dz; a=-0.5*dz; n=79; hd=(b-a)/(n+1); %********************************************************************* % THE IMPEDANCE MATRIX HAS A TOEPLITZ PROPERTY, THEREFORE ONLY % NM ELEMENTS NEED TO BE COMPUTED, AND THE MATRIX IS FILLED IN A % FORM THAT CAN BE SOLVED BY A TOEPLITZ MATRIX SOLVING %********************************************************************** for I = 1: nm zn=hl-(I-0.5)*dz; za1=zn-zm+a; %****************************************************************** % FAST ALGORITHM FORM OF THE SIMPSON'S INTEGRAL ROUTINE %****************************************************************** recgp=sqrt(ra*ra+za1*za1); cgp1=exp(-j*beta*recgp)*((1.0+j*beta*recgp)*(2.0*recgp*recgp-3.0*ra*ra)+(beta*ra*recgp)^2)/(2.0*beta*recgp^5); zb1=zn-zm+b; roc=sqrt(ra*ra+zb1*zb1); cgp2=exp(-j*beta*roc)*((1.0+j*beta*roc)*(2.0*roc*roc-3.0*ra*ra)+(beta*ra*roc)^2)/(2.0*beta*roc^5); crt=cgp1+cgp2; for k = 1: n xk=a+k*hd; zx1=zn-zm+xk; r=sqrt(ra*ra+zx1*zx1); cgp3=exp(-j*beta*r)*((1.0+j*beta*r)*(2.0*r*r-3.0*ra*ra)+(beta*ra*r)^2)/(2.0*beta*r^5); if mod(k,2)~=0 crt=crt+4.0*cgp3; else crt=crt+2.0*cgp3; end end crt=crt*hd*0.33333; zmn(I)=crt; if I~=1 zmn(nm+I-1)=crt; end end rb=2.3*ra; tlab=2.0*log(2.3); for i = 1: nm zi=hl-(i-0.5)*dz; r1=beta*sqrt(zi*zi+ra*ra); r2=beta*sqrt(zi*zi+rb*rb); j=0.0+j*1.0; if iex==1 cga(i)=-j*beta^2/(eta*tlab)*(expm(-j*r1)/r1-expm(-j*r2)/r2); else if i~=nmh cga(i)=0; else cga(i)=-j*beta/(eta*dz); end end end lolo = cga; %****************************************************************** % INPUT: % (C)A(2*M - 1) THE FIRST ROW OF THE T-MATRIX FOLLOWED BY % ITS FIRST COLUMN BEGINNING WITH THE SECOND % ELEMENT. ON RETURN A IS UNALTERED. % (C)B(M) THE RIGHT HAND SIDE VECTOR B. % (C)WA(2*M-2) A WORK AREA VECTOR % (I)M ORDER OF MATRIX A. % OUTPUT: % (C)B(M) THE SOLUTION VECTOR. % PURPOSE: % SOLVE A SYSTEM OF EQUATIONS DESCRIBED BY A TOEPLITZ MATRIX. % A * X = B % ****************************************************************** t2=length(zmn); t1=length(wa); m=nm; a=zmn; b=cga; wa=wa; m=nm; a1=a; a2=a(m+1:t2); t=b; c1=wa; c2=wa(m-1:t1); r1=a1(1); t(1)=b(1)/r1; if m ~=1 for n=2:m n1=n-1; n2=n-2; r5=a2(n1); r6=a1(n); if n ~= 2 c1(n1)=r2; for i1=1:n2 i2=n-i1; r5=r5+a2(i1)*c1(i2); r6=r6+a1(i1+1)*c2(i1); end end r2=-r5/r1; r3=-r6/r1; r1=r1+r5*r3; if n~=2 r6=c2(1); c2(n1)=0.0+j*0.0; for i1=2:n1 r5=c2(i1); c2(i1)=c1(i1)*r3+r6; c1(i1)=c1(i1)+r6*r2; r6=r5; end end c2(1)=r3; r5=0.0+j*0.0; for i1=1:n1 i2=n-i1; r5=r5+a2(i1)*t(i2); end r6=(b(n)-r5)/r1;