% Vivaldi_OpenEMS_Start.m % % Starter model for a CPW-fed, single-copper-layer, planar broadband % Vivaldi antenna on FR4, for the Studienprojekt "Design einer % Breitband-Vivaldi-Antenne" (Prof. Passoke). % % This script reuses the overall structure of the supplied Dipol_Passoke.m % (setup / substrate / metal / port / mesh / NF2FF / run / evaluate) and % replaces the straight dipole arms with an exponentially tapered Vivaldi % flare, built with AddPolygon (the P2 example that was commented out in % Dipol_Passoke.m is the template this is based on). % % STATUS: first-guess starting geometry, NOT yet optimized. Dimensions % come from the "Vivaldi antenna initial design" section of the project % strategy document. You must run this, check S11/Zin, and then sweep % parameters (Section 6 of the strategy doc) to converge on a design with % S11 < -10 dB from 850 MHz upward. % % ASSUMPTIONS made here that you should confirm with Prof. Passoke: % - Design/simulation band assumed 850 MHz - 3 GHz (upper limit not % given in the task sheet). % - Max board size assumed to fit within ~200 mm x 150 mm. % - CPW-fed single-layer design chosen instead of a double-sided % microstrip-to-slotline balun (task says a balun "may" be needed, % not that it is mandatory) - confirm this choice is acceptable. close all clear clc %% ------------------------------------------------------------------ %% 1) Constants & frequency setup %% ------------------------------------------------------------------ physical_constants; unit = 1e-3; % all lengths below are in mm f_low = 850e6; % required lower band edge f_high = 2000e6; % ASSUMED upper band edge -> confirm with professor f0 = (f_low + f_high)/2; % Gaussian excitation center frequency fc = (f_high - f_low)/2 + 0.3e9; % 20 dB corner freq, with margin below f_low % so the excitation spectrum still has % useful energy right down at 850 MHz FDTD = InitFDTD('NrTs', 3500000); % INCREASED from 2000000: the last run % (taper.L=160, ya=70, R=0.014) hit this % cap before reaching -50dB decay. % timesteps than a simple dipole to % decay fully - watch the end criteria % in the solver log; increase further % if it stops before decaying enough FDTD = SetGaussExcite(FDTD, f0, fc); BC = {'MUR' 'MUR' 'MUR' 'MUR' 'MUR' 'MUR'}; % start with MUR like the % dipole example; switch to % PML_8 if you see far-field % artifacts from reflections FDTD = SetBoundaryCond(FDTD, BC); %% ------------------------------------------------------------------ %% 2) Substrate parameters %% ------------------------------------------------------------------ % Air (kappa computed at design center frequency, same pattern as the % dipole starter file) luft.epsR = 1; luft.kappa = 0; % lossless air % FR4 substrate: epsR and tan delta from the datasheet / task sheet. % kappa = tan(delta) * 2*pi*f*eps0*epsR (constant-kappa approximation, % evaluated at the band center f0 - valid over a moderate bandwidth; % mention this approximation explicitly in your report) substrate.epsR = 4.3; substrate.tand = 0.023; substrate.kappa = substrate.tand * 2*pi*f0 * substrate.epsR * EPS0; substrate.thickness = 1; % mm, given by the task/datasheet %% ------------------------------------------------------------------ %% 3) Vivaldi geometry parameters (first-guess values, see strategy doc) %% ------------------------------------------------------------------ board.length = 275; % mm, taper direction (y-axis here) -- fits taper.L=240 % (throat_run 15 + 240 = 255, +20 mm margin). board.width = 160; % mm, transverse direction (x-axis here) -- fits % taper.ya=60 (2*(60+20) mm margin). feed.throat_run = 15; % mm, straight CPW section before the taper starts feed.R_ref = 50; % ohm, target system/port impedance (SMA/CPW) % CPW dimensions defined BEFORE the taper block on purpose: taper.y0 below % is derived from these so the ground plane's inner edge is geometrically % continuous at the throat/taper junction (see bug note below). cpw.center_w = 3.0; % mm, CPW center conductor width (feed section) cpw.gap = 0.4; % mm, CPW gap, center conductor to ground % NOTE: cpw.center_w / cpw.gap are analytic-formula starting guesses for % ~50 ohm on 1 mm FR4 (er=4.3). Finite ground width and the taper make % analytic CPW formulas inaccurate here - verify/tune these two numbers % against the simulated Zin, don't trust them blindly. taper.L = 240; % mm, taper length (throat -> aperture) -- BEST result % from the automated sweep (Vivaldi_Sweep.m): 6 combos % of taper.L in {160,200,240} x taper.ya in {50,60} were % tested. L=240/ya=60 had the best (least bad) worst- % case S11 across 850-2000 MHz (-3.18 dB), though NONE % of the 6 combos reached the -10 dB target continuously. % This is the production run to confirm that sweep % result with a finer mesh and full NrTs. % BUG FIX (see Report_Figures / KI-Nutzung notes): taper.y0 was previously % hardcoded to 1.0 mm, independent of cpw.center_w/cpw.gap. Since the feed % section's ground plane starts at x = cpw.center_w/2 + cpw.gap = 1.9 mm, % but the taper's ground inner edge started at x = 1.0 mm at the very same % y-position (y = feed.throat_run), the ground copper jumped inward past % the end of the 1.5 mm-half-width center conductor at that seam - i.e. % the ground plane physically overlapped/shorted the center conductor's % open end right at the throat/taper junction. This is almost certainly % the root cause of the S11 staying far from -10 dB and of the simulation % instabilities (oscillating/growing energy) seen when other dimensions % were changed, since changing board.length/taper.L moves the mesh lines % around this short and changes how it resonates. % FIX: start the taper's slot exactly where the feed section's gap ends, % so the ground plane is continuous (no overlap, no unintended gap). taper.y0 = cpw.center_w/2 + cpw.gap; % mm, HALF throat gap width at x=0 % (start of taper) = 1.9 mm taper.ya = 60.0; % mm, HALF aperture width at x=L (end of taper) -- % see taper.L note above (best sweep result). taper.R = 0.023; % 1/mm, taper rate - KEY OPTIMIZATION PARAMETER -- % REVERTED to the confirmed-stable value. Not touching % this again; taper.L is the lever being pushed now. %% ------------------------------------------------------------------ %% 4) Metal polygons: exponential taper curve generation %% ------------------------------------------------------------------ % y(x) = C1*exp(R*x) + C2 , 0 <= x <= L % with y(0) = y0 and y(L) = ya (see strategy doc Section 5 for derivation) N = 10; % kept at 10 (proven stable at L=240/ya=60 in the sweep) - not % risking a SmoothMesh hang on the production run by raising this. x_pts = linspace(0, taper.L, N); RL = taper.R * taper.L; C1 = (taper.ya - taper.y0) / (exp(RL) - 1); C2 = (taper.y0*exp(RL) - taper.ya) / (exp(RL) - 1); y_pts = C1*exp(taper.R*x_pts) + C2; % Place the throat at y-position 'feed.throat_run' measured from the % board's feed edge (y=0), so the full antenna occupies % y = 0 ... feed.throat_run + taper.L along the board. y_taper_abs = feed.throat_run + x_pts; % absolute y-coordinate along board % Build the closed polygon for ONE ground wing (the +x side): the ground % plane fills the board from the taper edge out to the board's side edge, % from the throat all the way to the aperture end. Point lists are % assembled as [x-row ; y-row], matching the AddPolygon(...,'CoordSystem',0) % convention used in the commented-out P2 example of Dipol_Passoke.m. % Inner boundary (the slot edge), going from throat to aperture: inner_x = y_pts; % x-coordinate = taper half-width inner_y = y_taper_abs; % y-coordinate = position along board % Outer boundary (board edge), going back from aperture to throat: outer_x = ones(1,N) * (board.width/2); outer_y = fliplr(y_taper_abs); Ppos = [ [inner_x, outer_x] ; [inner_y, outer_y] ]; % ground wing, +x side Pneg = [ -[inner_x, outer_x] ; [inner_y, outer_y] ]; % ground wing, -x side (mirrored) %% ------------------------------------------------------------------ %% 5) CSXCAD geometry %% ------------------------------------------------------------------ CSX = InitCSX(); % mesh air-box: must enclose the (larger) Vivaldi aperture + radiating % space; recompute/enlarge if you change board size SimBox = [board.width+300, board.length+300, 300]; mesh.x = [-SimBox(1)/2 SimBox(1)/2]; mesh.y = [-50 SimBox(2)-50]; mesh.z = [-SimBox(3)/2 SimBox(3)/2]; % substrate CSX = AddMaterial(CSX, 'substrate'); CSX = SetMaterialProperty(CSX, 'substrate', 'Epsilon', substrate.epsR, 'Kappa', substrate.kappa); start = [-board.width/2, 0, 0]; stop = [ board.width/2, board.length, substrate.thickness]; CSX = AddBox(CSX, 'substrate', 0, start, stop); % copper: CPW feed line (straight section from y=0 to throat) CSX = AddMetal(CSX, 'copper'); start = [-cpw.center_w/2, 0, substrate.thickness]; stop = [ cpw.center_w/2, feed.throat_run, substrate.thickness]; CSX = AddBox(CSX, 'copper', 10, start, stop); % copper: ground wings (tapered), both sides, built from the polygons above CSX = AddPolygon(CSX, 'copper', 10, 2, substrate.thickness, Ppos, 'CoordSystem', 0); CSX = AddPolygon(CSX, 'copper', 10, 2, substrate.thickness, Pneg, 'CoordSystem', 0); % also fill the ground on either side of the straight feed section % (from the CPW gap out to the board edge), so the CPW line has ground % right from the feed edge, not just starting at the throat Pfeed_pos = [ [cpw.center_w/2+cpw.gap, board.width/2, board.width/2, cpw.center_w/2+cpw.gap] ; ... [0, 0, feed.throat_run, feed.throat_run] ]; Pfeed_neg = [ -[cpw.center_w/2+cpw.gap, board.width/2, board.width/2, cpw.center_w/2+cpw.gap] ; ... [0, 0, feed.throat_run, feed.throat_run] ]; CSX = AddPolygon(CSX, 'copper', 10, 2, substrate.thickness, Pfeed_pos, 'CoordSystem', 0); CSX = AddPolygon(CSX, 'copper', 10, 2, substrate.thickness, Pfeed_neg, 'CoordSystem', 0); %% ------------------------------------------------------------------ %% 6) Feed / port %% ------------------------------------------------------------------ % Lumped port bridging the CPW gap at the feed edge (y=0), same mechanism % as the gap-bridging port in Dipol_Passoke.m, just oriented across the % CPW gap instead of along the dipole arms. start = [ cpw.center_w/2, 0, substrate.thickness]; stop = [ cpw.center_w/2+cpw.gap, 0.1, substrate.thickness]; % bridging the CPW gap [CSX, port] = AddLumpedPort(CSX, 15, 1, feed.R_ref, start, stop, [1 0 0], true); %% ------------------------------------------------------------------ %% 6b) Geometry sanity check (native Octave plot - no AppCSXCAD needed) %% ------------------------------------------------------------------ % AppCSXCAD/QCSXCAD is not installed on this Mac, so CSXGeomPlot() cannot % be used. This 2D plot draws the same copper polygons directly from the % arrays already computed above (feed line box, ground taper wings, feed % ground fill, and the port location) so you can catch overlap/gap bugs % BEFORE spending 10-20+ minutes on a full FDTD run. ALWAYS look at this % plot after changing any geometry parameter. figure('Name', 'Geometry check (copper layer, top view)'); hold on % feed line center conductor fill([-cpw.center_w/2 cpw.center_w/2 cpw.center_w/2 -cpw.center_w/2], ... [0 0 feed.throat_run feed.throat_run], [0.85 0.5 0.2]); % ground: feed-section fill + taper wings, both sides fill(Pfeed_pos(1,:), Pfeed_pos(2,:), [0.7 0.7 0.7]); fill(Pfeed_neg(1,:), Pfeed_neg(2,:), [0.7 0.7 0.7]); fill(Ppos(1,:), Ppos(2,:), [0.7 0.7 0.7]); fill(Pneg(1,:), Pneg(2,:), [0.7 0.7 0.7]); % port location (bridges the CPW gap at y=0) plot([cpw.center_w/2 cpw.center_w/2+cpw.gap], [0.05 0.05], 'r-', 'LineWidth', 4); % board outline plot([-board.width/2 board.width/2 board.width/2 -board.width/2 -board.width/2], ... [0 0 board.length board.length 0], 'k--'); axis equal xlabel('x / mm'); ylabel('y / mm'); title('Geometry check: orange=feed conductor, grey=ground, red=port - look for overlaps/gaps'); grid on drawnow disp('>>> CHECK THE GEOMETRY FIGURE NOW. Only continue if the orange feed'); disp(' line and the grey ground planes do NOT overlap, and the red port'); disp(' bar sits cleanly inside the gap between them.'); % NOTE: this check now runs BEFORE mesh generation (moved here on purpose) % so you see it quickly, without waiting through SmoothMesh, which can be % slow/hang for larger geometries. %% ------------------------------------------------------------------ %% 7) Mesh %% ------------------------------------------------------------------ % EXPLORATION-PHASE MESH: deliberately coarser than a final production % mesh (divisors reduced from /50,/20 to /20,/10) per the project's own % Schnellplan ("use a coarser mesh for exploratory sweep runs"). Larger % board/aperture sizes were taking 25+ minutes without fully converging % at the finer resolution - this trades some accuracy for much faster % iteration while we search for the right taper.L/ya/R combination. % IMPORTANT: once dimensions are chosen, switch these back to /50 and /20 % (or finer) and re-run once for the final, report-quality result. mesh = DetectEdges(CSX, mesh, 'ExcludeProperty', {'copper'}); mesh = DetectEdges(CSX, mesh, 'SetProperty', {'copper'}, '2D_Metal_Edge_Res', c0/(f_high)/unit/35); % resolution raised from /20 to /35 for this production/confirmation run % (finer than the exploration sweep, still short of the original /50 to % avoid re-triggering a SmoothMesh hang at this larger board size). % ROBUSTNESS FIX: pre-thin the raw DetectEdges line set BEFORE calling % SmoothMesh. SmoothMesh's recursive line-smoothing search has repeatedly % hung (100% CPU, 10-30+ min, zero output) whenever the taper curve's N % sample points happen to land close together in absolute position - % this has now happened with several different, unrelated parameter % changes (taper.L, board.width, taper.R), so it is SmoothMesh's own % algorithm that is fragile against dense/close line clusters, not any % one specific bug in our geometry. Thinning the raw line set first (any % two lines closer than pre_tol mm are merged) removes those clusters % before the expensive recursive step ever sees them. pre_tol = 0.15; % mm (still well below the smallest real feature, the 0.4 mm CPW gap) mesh.x = sort(mesh.x); mesh.x = mesh.x([true diff(mesh.x) > pre_tol]); mesh.y = sort(mesh.y); mesh.y = mesh.y([true diff(mesh.y) > pre_tol]); mesh.z = sort(mesh.z); mesh.z = mesh.z([true diff(mesh.z) > pre_tol]); mesh = SmoothMesh(mesh, c0/(f_high)/unit/15); % raised from /10 to /15 % BUG FIX: floating-point rounding in the exponential taper formula can % produce two mesh lines that are meant to be identical (e.g. both at % x = cpw.center_w/2+cpw.gap) but differ by ~1e-16 mm due to roundoff. % openEMS then builds a real (if absurdly tiny) cell between them, which % collapses the FDTD timestep to ~1e-20 s and makes the run meaningless. % Fix: snap/merge any mesh lines closer than 1e-6 mm (still far finer than % any real feature) on each axis before handing the mesh to openEMS. merge_tol = 1e-2; % mm (10 micron; still 40x below the smallest real feature, the 0.4 mm gap) mesh.x = sort(mesh.x); mesh.x = mesh.x([true diff(mesh.x) > merge_tol]); mesh.y = sort(mesh.y); mesh.y = mesh.y([true diff(mesh.y) > merge_tol]); mesh.z = sort(mesh.z); mesh.z = mesh.z([true diff(mesh.z) > merge_tol]); CSX = DefineRectGrid(CSX, unit, mesh); %% ------------------------------------------------------------------ %% 8) NF2FF box %% ------------------------------------------------------------------ start = [mesh.x(2) mesh.y(2) mesh.z(2)]; stop = [mesh.x(end-3) mesh.y(end-3) mesh.z(end-3)]; [CSX, nf2ff] = CreateNF2FFBox(CSX, 'nf2ff', start, stop); %% ------------------------------------------------------------------ %% 9) Run simulation %% ------------------------------------------------------------------ Sim_Path = 'Vivaldi01'; Sim_CSX = 'VIV01.xml'; [status, message, messageid] = rmdir(Sim_Path, 's'); [status, message, messageid] = mkdir(Sim_Path); WriteOpenEMS([Sim_Path '/' Sim_CSX], FDTD, CSX); % CSXGeomPlot([Sim_Path '/' Sim_CSX]); % ALWAYS check this plot before running - % confirm the taper/ground/feed look % correct before spending simulation time openEMS_opts = '--debug-PEC'; RunOpenEMS(Sim_Path, Sim_CSX, openEMS_opts); %% ------------------------------------------------------------------ %% 10) Post-processing: S11 / Zin %% ------------------------------------------------------------------ freq = linspace(f_low*0.7, f_high, 601); % start below f_low to see the % full roll-off of the match port = calcPort(port, Sim_Path, freq); Zin = port.uf.tot ./ port.if.tot; s11 = port.uf.ref ./ port.uf.inc; figure plot(freq/1e6, real(Zin), 'k-', 'Linewidth', 2); hold on; grid on plot(freq/1e6, imag(Zin), 'r--', 'Linewidth', 2); title('Feed point impedance'); xlabel('frequency f / MHz'); ylabel('impedance Z_{in} / Ohm'); legend('real', 'imag'); figure plot(freq/1e6, 20*log10(abs(s11)), 'k-', 'Linewidth', 2); grid on hold on plot([f_low f_low]/1e6, ylim, 'b--'); % marker at the required 850 MHz edge yline(-10, 'g--'); % -10 dB target line title('Reflection coefficient S_{11}'); xlabel('frequency f / MHz'); ylabel('|S_{11}| / dB'); drawnow %% ------------------------------------------------------------------ %% 11) Post-processing: far-field (NF2FF), directivity, efficiency %% ------------------------------------------------------------------ f_res_ind = find(20*log10(abs(s11)) == min(20*log10(abs(s11)))); f_res = freq(f_res_ind); disp('calculating far field at phi=[0 90] deg...'); nf2ff = CalcNF2FF(nf2ff, Sim_Path, f_res, [-180:2:180]*pi/180, [0 90]*pi/180); disp(['radiated power: Prad = ' num2str(nf2ff.Prad) ' Watt']); disp(['directivity: Dmax = ' num2str(nf2ff.Dmax) ' (' num2str(10*log10(nf2ff.Dmax)) ' dBi)']); disp(['efficiency: nu_rad = ' num2str(100*nf2ff.Prad./port.P_inc(f_res_ind)) ' %']); figure polarFF(nf2ff, 'xaxis', 'theta', 'param', [1 2], 'normalize', 1); figure plotFFdB(nf2ff, 'xaxis', 'theta', 'param', [1 2]); drawnow disp('calculating 3D far field pattern and dumping to vtk (use Paraview to visualize)...'); thetaRange = (0:2:180); phiRange = (0:2:360) - 180; nf2ff = CalcNF2FF(nf2ff, Sim_Path, f_res, thetaRange*pi/180, phiRange*pi/180, 'Verbose', 1, 'Outfile', '3D_Pattern.h5'); figure plotFF3D(nf2ff, 'logscale', -20); E_far_normalized = nf2ff.E_norm{1} / max(nf2ff.E_norm{1}(:)) * nf2ff.Dmax; DumpFF2VTK([Sim_Path '/3D_Pattern.vtk'], E_far_normalized, thetaRange, phiRange, 'scale', 1e-3); %% ------------------------------------------------------------------ %% NEXT STEPS (see strategy doc Section 6 for full detail) %% ------------------------------------------------------------------ % 1. Check CSXGeomPlot output BEFORE trusting the simulation - confirm % the taper shape, ground wings, and feed gap look geometrically % sensible (no overlaps, no gaps you didn't intend). % 2. Look at the S11 plot: find the lowest frequency where it crosses % below -10 dB. If that's above 850 MHz, increase taper.L and/or % taper.ya first (biggest levers), then re-run. % 3. Once S11 looks reasonable, verify Zin real part is close to 50 ohm % and imaginary part is small, over your target band. % 4. Sweep taper.R, cpw.center_w, cpw.gap, feed.throat_run one at a time, % tracking how the -10 dB crossing frequency and Zin move. % 5. After the design converges, refine the mesh (smaller % 2D_Metal_Edge_Res / SmoothMesh target) and re-run once to confirm % the result is mesh-converged, not a meshing artifact.