% % bi-quad antenna % 2.6GHz % % openEMS 64bit -- version v0.0.35-64-g3d684ba % octave -- version 5.2.0 % (C) 2020 Yellow Rabbit close all clear clc output_precision(5); threads_num = 2; %debug = ' --no-simulation'; %debug = ' --debug-PEC --no-simulation'; debug = ' '; show = 1; %% prepare simulation folder Sim_Path = tilde_expand('~/tmp/b26'); Sim_CSX = 'bi_quad_ant.xml'; confirm_recursive_rmdir(0); [status, message, messageid] = rmdir( Sim_Path, 's' ); % clear previous directory [status, message, messageid] = mkdir( Sim_Path ); % create empty simulation folder %% setup the simulation physical_constants; unit = 1e-3; % all length in mm % frequency range of interest f_start = 2.2e9; f0 = 2.6e9; f_stop = 2*f0 - f_start; freq = linspace(f_start, f_stop, 401); end_crit = 1e-6; lambda = c0/f0/unit; % wire wh = 2.2; ww = 1.3; Rc = 8/2; Ri = Rc - ww/2; Ro = Rc + ww/2; Rc_proj = Rc/sqrt(2); port_length = 1 + ww; port_proj = port_length * sqrt(2); dR = Rc*sqrt(2) - Rc - port_length/2; quad_size = lambda/4 * 0.973; delta_port = (2*pi*Rc)/4; quad_small = (2*quad_size - dR*sqrt(2) - 2*delta_port)/2; quad_big = (2*quad_size + dR*sqrt(2) - 2*delta_port)/2; reflector_off = lambda/8.25 + wh/2; reflector_thickness = wh; reflector_size = quad_size*4 * 1.1; quad_mesh = ww/2; port_mesh = port_length/6; max_res = lambda/15; reflector_mesh = lambda/20; Feed_R = 50; % size of the simulation box SimBox = [300 800 800]; %% setup FDTD parameter & excitation function FDTD = InitFDTD('endCriteria', end_crit); FDTD = SetGaussExcite(FDTD, 0.5*(f_start + f_stop), 0.5*(f_stop - f_start)); BC = {'PML_8' 'PML_8' 'PML_8' 'PML_8' 'PML_8' 'PML_8'}; % boundary conditions FDTD = SetBoundaryCond(FDTD, BC); %% setup CSXCAD geometry & mesh CSX = InitCSX(); % create fixed lines for the antenna port mesh.x = [-wh -wh/2 0 wh/2 wh]; mesh.y = [-Rc_proj-ww/2 -Rc_proj+ww/2 0 Rc_proj-ww/2 Rc_proj+ww/2]; mesh.z = [-port_length/2-ww/2 -port_length/2+ww/2 0 port_length/2-ww/2 port_length/2+ww/2]; mesh = SmoothMesh(mesh, port_mesh, 1.3); % create fixed lines for the antenna outline mesh.x = [mesh.x -reflector_off-reflector_thickness -reflector_off]; mesh.y = [mesh.y -(3*Rc + quad_big + quad_small)/sqrt(2)+Rc_proj-Rc-ww/2 ... -(3*Rc + quad_big + quad_small)/sqrt(2)+Rc_proj-Rc+ww/2 ... (3*Rc + quad_big + quad_small)/sqrt(2)-Rc_proj+Rc-ww/2 ... (3*Rc + quad_big + quad_small)/sqrt(2)-Rc_proj+Rc+ww/2]; mesh.z = [mesh.z -(Rc + quad_small)/sqrt(2)+Rc_proj-Rc-ww/2 ... -(Rc + quad_small)/sqrt(2)+Rc_proj-Rc+ww/2 ... (Rc + quad_small)/sqrt(2)-Rc_proj+Rc-ww/2 ... (Rc + quad_small)/sqrt(2)-Rc_proj+Rc+ww/2]; mesh = SmoothMesh(mesh, quad_mesh, 1.3); % create reflector mesh mesh.y = [mesh.y -reflector_size/2 reflector_size/2]; mesh.z = [mesh.z -reflector_size/2 reflector_size/2]; mesh = SmoothMesh(mesh, reflector_mesh, 1.3); % add air box mesh.x = [mesh.x -SimBox(1)/4 SimBox(1)*3/4]; mesh.y = [mesh.y -SimBox(2)/2 SimBox(2)/2]; mesh.z = [mesh.z -SimBox(3)/2 SimBox(3)/2]; mesh = SmoothMesh(mesh, max_res, 1.5); CSX = DefineRectGrid(CSX, unit, mesh); % create a biquad CSX = AddMaterial(CSX,'metal'); #CSX = AddMetal(CSX,'metal'); %% Copper CSX = SetMaterialProperty( CSX, 'metal', 'Kappa', 56e6); % inner CSX = AddBox(CSX, 'metal', 10, [-wh/2 -ww/2 0], [wh/2 ww/2 quad_big], 'Transform', {'Rotate_X', -pi/4, 'Translate', ... ['0,' num2str(Rc_proj) ',' num2str(Rc_proj - dR)]}); CSX = AddBox(CSX, 'metal', 10, [-wh/2 -ww/2 0], [wh/2 ww/2 -quad_big], 'Transform', {'Rotate_X', pi/4, 'Translate', ... ['0,' num2str(Rc_proj) ',' num2str(-Rc_proj + dR)]}); CSX = AddBox(CSX, 'metal', 10, [-wh/2 -ww/2 0], [wh/2 ww/2 quad_big], 'Transform', {'Rotate_X', pi/4, 'Translate', ... ['0,' num2str(-Rc_proj) ',' num2str(Rc_proj -dR)]}); CSX = AddBox(CSX, 'metal', 10, [-wh/2 -ww/2 0], [wh/2 ww/2 -quad_big], 'Transform', {'Rotate_X', -pi/4, 'Translate', ... ['0,' num2str(-Rc_proj) ',' num2str(-Rc_proj + dR)]}); % outer CSX = AddBox(CSX, 'metal', 10, [-wh/2 -ww/2 0], [wh/2 ww/2 -quad_small], 'Transform', {'Rotate_X', pi/4, 'Translate', ... ['0,' num2str(3*Rc_proj + quad_big/sqrt(2)) ',' num2str((quad_small + Rc)/sqrt(2))]}); CSX = AddBox(CSX, 'metal', 10, [-wh/2 -ww/2 0], [wh/2 ww/2 -quad_small], 'Transform', {'Rotate_X', -pi/4, 'Translate', ... ['0,' num2str(-3*Rc_proj - quad_big/sqrt(2)) ',' num2str((quad_small + Rc)/sqrt(2))]}); CSX = AddBox(CSX, 'metal', 10, [-wh/2 -ww/2 0], [wh/2 ww/2 quad_small], 'Transform', {'Rotate_X', pi/4, 'Translate', ... ['0,' num2str(-3*Rc_proj - quad_big/sqrt(2)) ',' num2str(-(quad_small + Rc)/sqrt(2))]}); CSX = AddBox(CSX, 'metal', 10, [-wh/2 -ww/2 0], [wh/2 ww/2 quad_small], 'Transform', {'Rotate_X', -pi/4, 'Translate', ... ['0,' num2str(3*Rc_proj + quad_big/sqrt(2)) ',' num2str(-(quad_small + Rc)/sqrt(2))]}); % curve corners p(1, 1) = Ri; p(2, 1) = -wh/2; p(1, 2) = Ri; p(2, 2) = wh/2; p(1, 3) = Ro; p(2, 3) = wh/2; p(1, 4) = Ro; p(2, 4) = -wh/2; CSX = AddRotPoly(CSX, 'metal', 10, 'y', p, 'x', [3*pi/4 5*pi/4], 'Transform', ... {'Translate', ['0,0,' num2str(Rc + port_length/2)]}); CSX = AddRotPoly(CSX, 'metal', 10, 'y', p, 'x', [7*pi/4 9*pi/4], 'Transform', ... {'Translate', ['0,0,' num2str(-Rc - port_length/2)]}); CSX = AddRotPoly(CSX, 'metal', 10, 'y', p, 'x', [5*pi/4 7*pi/4], 'Transform', ... {'Translate', ['0,' num2str((quad_big + quad_small + 2*Rc)/sqrt(2)) ',0']}); CSX = AddRotPoly(CSX, 'metal', 10, 'y', p, 'x', [pi/4 3*pi/4], 'Transform', ... {'Translate', ['0,' num2str(-(quad_big + quad_small + 2*Rc)/sqrt(2)) ',0']}); CSX = AddRotPoly(CSX, 'metal', 10, 'y', p, 'x', [7*pi/4 9*pi/4], 'Transform', ... {'Translate', ['0,' num2str(-(quad_big + 2*Rc)/sqrt(2)) ',' ... num2str((Rc + quad_small)/sqrt(2) - Rc_proj)]}); CSX = AddRotPoly(CSX, 'metal', 10, 'y', p, 'x', [7*pi/4 9*pi/4], 'Transform', ... {'Translate', ['0,' num2str((quad_big + 2*Rc)/sqrt(2)) ',' ... num2str((Rc + quad_small)/sqrt(2) - Rc_proj)]}); CSX = AddRotPoly(CSX, 'metal', 10, 'y', p, 'x', [3*pi/4 5*pi/4], 'Transform', ... {'Translate', ['0,' num2str(-(quad_big + 2*Rc)/sqrt(2)) ',' ... num2str(-(Rc + quad_small)/sqrt(2) + Rc_proj)]}); CSX = AddRotPoly(CSX, 'metal', 10, 'y', p, 'x', [3*pi/4 5*pi/4], 'Transform', ... {'Translate', ['0,' num2str((quad_big + 2*Rc)/sqrt(2)) ',' ... num2str(-(Rc + quad_small)/sqrt(2) + Rc_proj)]}); % reflector plane CSX = AddMetal(CSX,'reflector'); start = [-reflector_off-reflector_thickness -reflector_size/2 -reflector_size/2]; stop = [-reflector_off reflector_size/2 reflector_size/2]; CSX = AddBox(CSX, 'reflector', 10, start, stop); %% apply the excitation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% start = [0 0 -port_length/2+ww/2]; stop = [0 0 port_length/2-ww/2]; [CSX port] = AddLumpedPort(CSX, 10, 0, Feed_R, start, stop, [0 0 1], true); %% nf2ff calc start = [mesh.x(9) mesh.y(9) mesh.z(9)]; stop = [mesh.x(end-8) mesh.y(end-8) mesh.z(end-8)]; [CSX nf2ff] = CreateNF2FFBox(CSX, 'nf2ff', start, stop); %% write openEMS compatible xml-file WriteOpenEMS([Sim_Path '/' Sim_CSX], FDTD, CSX); %% show the structure if (show) CSXGeomPlot([Sim_Path '/' Sim_CSX], ['--export-polydata-vtk=' Sim_Path]); end %% run openEMS Settings.LogFile = [Sim_Path '/openEMS.log']; options = ['--engine=multithreaded --numThreads=' num2str(threads_num) debug]; RunOpenEMS(Sim_Path, Sim_CSX, options, Settings); % Postprocessing & do the plots freq = linspace(f_start, f_stop, 201); port = calcPort(port, Sim_Path, freq); Zin = port.uf.tot ./ port.if.tot; S11 = port.uf.ref ./ port.uf.inc; %% swr & resonant swr = (1 + abs(S11)) ./ (1 - abs(S11)); [swr_min idx] = min(swr); [Xmin idx1] = min(abs(imag(Zin))); swr_fr = freq(idx); res_fr = freq(idx1); Zr_swr = Zin(idx); Zr_res = Zin(idx1); disp('Minimum SWR frequency:'); disp(swr_fr); disp('Resonant frequency (jX=0):'); disp(res_fr); disp('Impedance at minimum SWR:'); disp(Zr_swr); disp('Impedance at resonant frequency:'); disp(Zr_res); disp('Minimum SWR:'); disp(swr_min); % S11 plot(freq/1e6, 20*log10(abs(S11))); grid on title('Reflection coefficient'); xlabel('frequency f / MHz'); ylabel('S11 / dB'); print([Sim_Path '/s11.png']); close; % plot Zin figure plot(freq/1e6, real(Zin), 'k-', 'Linewidth', 2); grid on title('feed point impedance'); xlabel('frequency f / MHz'); ylabel('resistance Z_{in} / Ohm'); print([Sim_Path '/Zin-real.png']); close; figure plot(freq/1e6, imag(Zin), 'r--', 'Linewidth', 2); grid on title('feed point impedance'); xlabel('frequency f / MHz'); ylabel('reactance Z_{in} / Ohm'); print([Sim_Path '/Zin-imag.png']); close; %% calculate 3D far field pattern phiRange = -180:2.5:180; thetaRange = 0:2.5:180; nf2ff = CalcNF2FF(nf2ff, Sim_Path, f0, thetaRange*pi/180, phiRange*pi/180); disp( ['directivity: Dmax = ' num2str(10*log10(nf2ff.Dmax)) ' dBi'] ); % plot far-field pattern with Matlab %figure %plotFF3D(nf2ff, 'logscale', -20) %% disp( 'Dumping far-field pattern to vtk (use Paraview to visualize)...' ); DumpFF2VTK([Sim_Path '/Bi_Quad_Pattern.vtk'], nf2ff.E_norm{1} / max(nf2ff.E_norm{1}(:)) * nf2ff.Dmax, thetaRange, phiRange, 'scale', 0.05); % vim: expandtab: sw=4 ts=4: