1# -*- coding: utf-8 -*- 2""" 3 Bent Patch Antenna Tutorial 4 5 Tested with 6 - python 3.4 7 - openEMS v0.0.33+ 8 9 (C) 2016 Thorsten Liebig <thorsten.liebig@gmx.de> 10 11""" 12 13### Import Libraries 14import os, tempfile 15from pylab import * 16from mpl_toolkits.mplot3d import Axes3D 17 18from CSXCAD import CSXCAD 19 20from openEMS.openEMS import openEMS 21from openEMS.physical_constants import * 22 23 24### Setup the simulation 25Sim_Path = os.path.join(tempfile.gettempdir(), 'Bent_Patch') 26 27post_proc_only = False 28 29unit = 1e-3 # all length in mm 30 31f0 = 2.4e9 # center frequency, frequency of interest! 32lambda0 = round(C0/f0/unit) # wavelength in mm 33fc = 0.5e9 # 20 dB corner frequency 34 35# patch width in alpha-direction 36patch_width = 32 # resonant length in alpha-direction 37patch_radius = 50 # radius 38patch_length = 40 # patch length in z-direction 39 40#substrate setup 41substrate_epsR = 3.38 42substrate_kappa = 1e-3 * 2*pi*2.45e9 * EPS0*substrate_epsR 43substrate_width = 80 44substrate_length = 90 45substrate_thickness = 1.524 46substrate_cells = 4 47 48#setup feeding 49feed_pos = -5.5 #feeding position in x-direction 50feed_width = 2 #feeding port width 51feed_R = 50 #feed resistance 52 53# size of the simulation box 54SimBox_rad = 2*100 55SimBox_height = 1.5*200 56 57### Setup FDTD parameter & excitation function 58FDTD = openEMS(CoordSystem=1, EndCriteria=1e-4) # init a cylindrical FDTD 59f0 = 2e9 # center frequency 60fc = 1e9 # 20 dB corner frequency 61FDTD.SetGaussExcite(f0, fc) 62FDTD.SetBoundaryCond(['MUR', 'MUR', 'MUR', 'MUR', 'MUR', 'MUR']) # boundary conditions 63 64### Setup the Geometry & Mesh 65# init a cylindrical mesh 66CSX = CSXCAD.ContinuousStructure(CoordSystem=1) 67FDTD.SetCSX(CSX) 68mesh = CSX.GetGrid() 69mesh.SetDeltaUnit(unit) 70 71### Setup the geometry using cylindrical coordinates 72# calculate some width as an angle in radiant 73patch_ang_width = patch_width/(patch_radius+substrate_thickness) 74substr_ang_width = substrate_width/patch_radius 75feed_angle = feed_pos/patch_radius 76 77# create patch 78patch = CSX.AddMetal('patch') # create a perfect electric conductor (PEC) 79start = [patch_radius+substrate_thickness, -patch_ang_width/2, -patch_length/2 ] 80stop = [patch_radius+substrate_thickness, patch_ang_width/2, patch_length/2 ] 81patch.AddBox(priority=10, start=start, stop=stop) # add a box-primitive to the metal property 'patch' 82FDTD.AddEdges2Grid(dirs='all', properties=patch) 83 84# create substrate 85substrate = CSX.AddMaterial('substrate', epsilon=substrate_epsR, kappa=substrate_kappa ) 86start = [patch_radius , -substr_ang_width/2, -substrate_length/2] 87stop = [patch_radius+substrate_thickness, substr_ang_width/2, substrate_length/2] 88substrate.AddBox(start=start, stop=stop) 89FDTD.AddEdges2Grid(dirs='all', properties=substrate) 90 91# save current density oon the patch 92jt_patch = CSX.AddDump('Jt_patch', dump_type=3, file_type=1) 93start = [patch_radius+substrate_thickness, -substr_ang_width/2, -substrate_length/2] 94stop = [patch_radius+substrate_thickness, +substr_ang_width/2, substrate_length/2] 95jt_patch.AddBox(start=start, stop=stop) 96 97# create ground 98gnd = CSX.AddMetal('gnd') # create a perfect electric conductor (PEC) 99start = [patch_radius, -substr_ang_width/2, -substrate_length/2] 100stop = [patch_radius, +substr_ang_width/2, +substrate_length/2] 101gnd.AddBox(priority=10, start=start, stop=stop) 102FDTD.AddEdges2Grid(dirs='all', properties=gnd) 103 104# apply the excitation & resist as a current source 105start = [patch_radius , feed_angle, 0] 106stop = [patch_radius+substrate_thickness, feed_angle, 0] 107port = FDTD.AddLumpedPort(1 ,feed_R, start, stop, 'r', 1.0, priority=50, edges2grid='all') 108 109### Finalize the Mesh 110# add the simulation domain size 111mesh.AddLine('r', patch_radius+np.array([-20, SimBox_rad])) 112mesh.AddLine('a', [-0.75*pi, 0.75*pi]) 113mesh.AddLine('z', [-SimBox_height/2, SimBox_height/2]) 114 115# add some lines for the substrate 116mesh.AddLine('r', patch_radius+np.linspace(0,substrate_thickness,substrate_cells)) 117 118# generate a smooth mesh with max. cell size: lambda_min / 20 119max_res = C0 / (f0+fc) / unit / 20 120max_ang = max_res/(SimBox_rad+patch_radius) # max res in radiant 121mesh.SmoothMeshLines(0, max_res, 1.4) 122mesh.SmoothMeshLines(1, max_ang, 1.4) 123mesh.SmoothMeshLines(2, max_res, 1.4) 124 125## Add the nf2ff recording box 126nf2ff = FDTD.CreateNF2FFBox() 127 128### Run the simulation 129if 0: # debugging only 130 CSX_file = os.path.join(Sim_Path, 'bent_patch.xml') 131 if not os.path.exists(Sim_Path): 132 os.mkdir(Sim_Path) 133 CSX.Write2XML(CSX_file) 134 os.system(r'AppCSXCAD "{}"'.format(CSX_file)) 135 136 137if not post_proc_only: 138 FDTD.Run(Sim_Path, verbose=3, cleanup=True) 139 140### Postprocessing & plotting 141f = np.linspace(max(1e9,f0-fc),f0+fc,401) 142port.CalcPort(Sim_Path, f) 143Zin = port.uf_tot / port.if_tot 144s11 = port.uf_ref/port.uf_inc 145s11_dB = 20.0*np.log10(np.abs(s11)) 146 147figure() 148plot(f/1e9, s11_dB) 149grid() 150ylabel('s11 (dB)') 151xlabel('frequency (GHz)') 152 153P_in = 0.5*np.real(port.uf_tot * np.conj(port.if_tot)) # antenna feed power 154 155# plot feed point impedance 156figure() 157plot( f/1e6, real(Zin), 'k-', linewidth=2, label=r'$\Re(Z_{in})$' ) 158grid() 159plot( f/1e6, imag(Zin), 'r--', linewidth=2, label=r'$\Im(Z_{in})$' ) 160title( 'feed point impedance' ) 161xlabel( 'frequency (MHz)' ) 162ylabel( 'impedance ($\Omega$)' ) 163legend( ) 164 165 166idx = np.where((s11_dB<-10) & (s11_dB==np.min(s11_dB)))[0] 167if not len(idx)==1: 168 print('No resonance frequency found for far-field calulation') 169else: 170 f_res = f[idx[0]] 171 theta = np.arange(-180.0, 180.0, 2.0) 172 print("Calculate NF2FF") 173 nf2ff_res_phi0 = nf2ff.CalcNF2FF(Sim_Path, f_res, theta, 0, center=np.array([patch_radius+substrate_thickness, 0, 0])*unit, read_cached=True, outfile='nf2ff_xz.h5') 174 175 figure(figsize=(15, 7)) 176 ax = subplot(121, polar=True) 177 E_norm = 20.0*np.log10(nf2ff_res_phi0.E_norm/np.max(nf2ff_res_phi0.E_norm)) + nf2ff_res_phi0.Dmax 178 ax.plot(np.deg2rad(theta), 10**(np.squeeze(E_norm)/20), linewidth=2, label='xz-plane') 179 ax.grid(True) 180 ax.set_xlabel('theta (deg)') 181 ax.set_theta_zero_location('N') 182 ax.set_theta_direction(-1) 183 ax.legend(loc=3) 184 185 phi = theta 186 nf2ff_res_theta90 = nf2ff.CalcNF2FF(Sim_Path, f_res, 90, phi, center=np.array([patch_radius+substrate_thickness, 0, 0])*unit, read_cached=True, outfile='nf2ff_xy.h5') 187 188 ax = subplot(122, polar=True) 189 E_norm = 20.0*np.log10(nf2ff_res_theta90.E_norm/np.max(nf2ff_res_theta90.E_norm)) + nf2ff_res_theta90.Dmax 190 ax.plot(np.deg2rad(phi), 10**(np.squeeze(E_norm)/20), linewidth=2, label='xy-plane') 191 ax.grid(True) 192 ax.set_xlabel('phi (deg)') 193 suptitle('Bent Patch Anteanna Pattern\nFrequency: {} GHz'.format(f_res/1e9), fontsize=14) 194 ax.legend(loc=3) 195 196 print( 'radiated power: Prad = {:.2e} Watt'.format(nf2ff_res_theta90.Prad[0])) 197 print( 'directivity: Dmax = {:.1f} ({:.1f} dBi)'.format(nf2ff_res_theta90.Dmax[0], 10*np.log10(nf2ff_res_theta90.Dmax[0]))) 198 print( 'efficiency: nu_rad = {:.1f} %'.format(100*nf2ff_res_theta90.Prad[0]/real(P_in[idx[0]]))) 199 200show() 201 202