minimal_dipole.py

import os
import math
import tempfile
import numpy as np
from CSXCAD import ContinuousStructure, AppCSXCAD_BIN
from CSXCAD.CSProperties import CSPropMetal
from openEMS import openEMS
from openEMS.physical_constants import C0

unit = 1e-3

# 100-500 MHz: f0-bw, f0+bw
f0 = 300e6
bw = 200e6

lambda_f_max = C0 / (f0 + bw) / unit
lambda_f_min = C0 / (f0 - bw) / unit

dipole_length = 285.8 / 2
dipole_gap = 1.0
dipole_wire_radius = 2.0

feed_resistance = 70.8
feed_radius = dipole_wire_radius / math.sqrt(2)

sim_box = np.array([1, 1, 1]) * 2.0 * lambda_f_min
max_res = math.floor(lambda_f_max / 20)

fdtd = openEMS(NrTS=100000, EndCriteria=1e-4)
fdtd.SetDiracExcite(lambda_f_max)
fdtd.SetBoundaryCond(["MUR", "MUR", "MUR", "MUR", "MUR", "MUR"])
csx = ContinuousStructure()
fdtd.SetCSX(csx)

mesh = csx.GetGrid()
mesh.SetDeltaUnit(unit)

mesh.AddLine("z", np.linspace(-dipole_gap / 2, dipole_gap / 2, 5))
mesh.AddLine("z",
             np.linspace(-dipole_length / 2 - 5 * dipole_wire_radius, -dipole_length / 2 + 5 * dipole_wire_radius, 11))
mesh.AddLine("z",
             np.linspace(dipole_length / 2 - 5 * dipole_wire_radius, dipole_length / 2 + 5 * dipole_wire_radius, 11))

mesh.AddLine("x", [-sim_box[0] / 2, 0, sim_box[0] / 2])
mesh.SmoothMeshLines("x", max_res, ratio=1.4)

mesh.AddLine("y", [-sim_box[1] / 2, 0, sim_box[1] / 2])
mesh.SmoothMeshLines("y", max_res, ratio=1.4)

mesh.AddLine("z", [-sim_box[2] / 2, 0, sim_box[2] / 2])
mesh.SmoothMeshLines("z", max_res, ratio=1.4)

arm1: CSPropMetal = csx.AddMetal("arm1")
arm1.AddWire([[0, 0], [0, 0], [-dipole_gap / 2, -dipole_length / 2]], radius=dipole_wire_radius)
arm1.SetColor("#ff0000", 50)

arm2: CSPropMetal = csx.AddMetal("arm2")
arm2.AddWire([[0, 0], [0, 0], [dipole_gap / 2, dipole_length / 2]], radius=dipole_wire_radius)
arm2.SetColor("#ff0000", 50)

feed = fdtd.AddLumpedPort(
    1,
    feed_resistance,
    [-feed_radius, -feed_radius, -dipole_gap / 2],
    [feed_radius, feed_radius, dipole_gap / 2],
    "z",
    1.0,
    priority=5,
)

output_dir = os.path.join(tempfile.gettempdir(), "dipole")
if not os.path.isdir(output_dir):
    os.mkdir(output_dir)
output_fn = os.path.join(output_dir, "dipole.xml")
csx.Write2XML(output_fn)
os.system('{} "{}"'.format(AppCSXCAD_BIN, output_fn))

fdtd.Run(output_dir, verbose=3, cleanup=True)

freq = np.linspace(f0 - bw, f0 + bw, 401)
feed.CalcPort(output_dir, freq)

Zin = feed.uf_tot / feed.if_tot
s11 = feed.uf_ref / feed.uf_inc
s11_dB = 20.0 * np.log10(np.abs(s11))

idx = np.where((s11_dB < -10.0) & (s11_dB == np.min(s11_dB)))[0]
if not len(idx) == 1:
    print("No resonance frequency found for far-field calculation!")
else:
    print(
        "Found resonance frequency at {} MHz with {} dB at {} Ohm".format(
            round(freq[idx][0] / 1e6, 1),
            round(s11_dB[idx][0], 1),
            round(np.real(Zin[idx])[0], 1),
        )
    )