RFSimulatorFEMWELL
The RFSimulatorFEMWELL is a FEM mode solver built around FEMWELL and skfem. With this tool we can calculate the RF mode supported by our structure. We can retrieve the propagation constant, and calculate the characteristic impedance, small signal S parameters and equivalent RLGC transmission line model parameters [[1], [3]].
- class imodulator.RFSimulator.RFSimulatorFEMWELL(device, simulation_window=None)
Base class for RF simulation of a PhotonicDevice. The role of the RFsimulator is to:
Warning
The algorithm to search for junctions is not active at the moment. It is a feature that will be implemented in the future.
generate the RF mesh;
apply RF simmetry planes;
Generate field visualizations of the RF modes;
Generate mesh visuals;
Solve for the RF modes;
Return small-signal S parameters;
Return RLGC parameters based on waveguide theory;
- __init__(device, simulation_window=None)
Initializes the simulator with the photonic device and the simulation window.
- Parameters:
device (
PhotonicDevice) – ThePhotonicDeviceobject to simulate.simulation_window (
Polygon|None(default:None)) – The simulation window to use for the RF simulation. If None, the entire device is simulated.
Note
The simulation window MUST be a rectangle. If not, the simulation will fail. Also beware that the definition of the simulation window is how you can make use of symmetry planes through the definition of metal boundaries or not. See compute_modes for more information on how to use the symmetry planes.
- compute_modes(frequency=10, voltage_idx=0, num_modes=1, order=1, metallic_boundaries=False, n_guess=4.0, return_modes=False, use_charge_transport_data=False)
Compute the electromagnetic RF modes of a photonic device at a given frequency. The modes are computed via FEMWELL.
- Parameters:
frequency (
float(default:10)) – The frequency at which to compute the modes. The frequency must be in GHz.voltage_idx (
int(default:0)) – The index of the voltage to use for the permittivity tensor.num_modes (
int(default:1)) – The number of modes to compute.order (
int(default:1)) – Order of the basis functions to use in the EM solver.metallic_boundaries (
list|str|bool(default:False)) – The boundaries to treat as metallic. If False, no boundaries are treated as metallic. If True, all boundaries are treated as metallic. If a list of strings, the boundaries with the given names are treated as metallic. At the moment, the simulation window is treated as a square, therefore, the metallic boundaries can beleft`, ``right,topandbottomboundaries.n_guess (
float(default:4.0)) – Initial guess for the effective index.return_modes (
bool(default:False)) – Whether to return the computed modes.use_charge_transport_data (
bool(default:False)) – Whether to use the charge transport data to compute the permittivity tensor. Doing so will yield a \(\sigma(x,y,V)\). Make sure your mesh is appropriate for it.
- Return type:
Modes- Returns:
The computed modes if return_modes is True.
- get_RLGC(mode, i0=1.0, v0=1.0)
This function calculates the RLGC parameters calculated from waveguide circuit theory for a TEM field.
The \(i_0\) and \(v_0\) values provided should be the values of current and voltage across the signal and ground tracks of the equivalent transmission line. This is not always straightforward, particularly if you’re dealing with multiconductor transmission lines. It is reccomended to use this formalism only when you have a good equivalent transmission line of the present mode.
Please note that every field returned by FEMWELL is power normalized by orthogonality relations. Therefore, if you’re using a symmetry plane, you have to take into account that this code will return a value taking into account that HALF the field has power 1. Please adjust the values to your specific case/symmetry. The reccomended approach is to calculate once with the full structure and then find the proper scaling factors.
The RLGC parameters are calculated according to Marks and Williams [1]:
\[\begin{split}C = \frac{1}{|v_0|^2} \left[ \int_S \Re\{\epsilon_{RF}\} |E_t|^2dS - \int_S \Re\{\mu_{RF}\}|H_z|^2 dS \right] \\ L = \frac{1}{|i_0|^2} \left[ \int_S \Re\{\mu_{RF}\} |H_t|^2dS - \int_S \Re\{\epsilon_{RF}\}|E_z|^2 dS \right] \\ G = \frac{\omega}{|v_0|^2} \left[ \int_S \Im\{\epsilon_{RF}\} |E_t|^2dS + \int_S \Im\{\mu_{RF}\}|H_z|^2 dS \right] \\ R = \frac{\omega}{|i_0|^2} \left[ \int_S \Im\{\mu_{RF}\} |H_t|^2dS + \int_S \Im\{\epsilon_{RF}\}|E_z|^2 dS \right]\end{split}\]- Parameters:
mode (
Mode) – The mode object containing the electric and magnetic field data to be used for the calculation.i0 (
complex(default:1.0)) – The current value to use for the calculation.v0 (
complex(default:1.0)) – The voltage value to use for the calculation.
- Returns:
The resistance per unit length in Ohms per micrometer. L (float): The inductance per unit length in Henrys per micrometer. G (float): The conductance per unit length in Siemens per micrometer. C (float): The capacitance per unit length in Farads per micrometer.
- Return type:
R (float)
- get_S(gamma, Z, ZL=50.0, ZS=50.0, L=1000.0)
We calculate the S parameters based on transmission line theory Rizzi [3], starting from the ABCD matrix and going to S parameters:
\[\begin{split}S_{11} = \frac{AZ_L + B - CZ_L*Z_S - DZ_S}{AZ_L + B + CZ_LZ_S + DZ_S} \\ S_{22} = \frac{2\sqrt{Z_S Z_L}}{AZ_L + B + CZ_LZ_S + DZ_S}\end{split}\]- Parameters:
gamma (
Quantity) – The propagation constant. Must be in um**-1Z (
Quantity) – The characteristic impedance. Must be in Ohm.ZL (
float(default:50.0)) – The load impedance. Defaults to 50 Ohm.ZS (
float(default:50.0)) – The source impedance. Defaults to 50 Ohm.L (
float(default:1000.0)) – The length of the transmission line. Defaults to 1000 um.
- Returns:
The S11 parameter. S12: The S12 parameter
- Return type:
S11
- get_currents(mode)
This function calculates the currents in every polygon marked as
calculate_current = Truevia:\[i_0 = \int_{\partial \Omega} \mathbf{n} \times \mathbf{H} \cdot d\mathbf{s}\]Note that every field returned by FEMWELL is power normalized by orthogonality relations. Therefore, if you’re using a symmetry plane, you have to take into account that this code will return a value taking into account that HALF the field has power 1. Please adjust the values to your specific case/symmetry. The reccomended approach is to calculate once with the full structure and then find the proper scaling factors.
- Parameters:
mode (
Mode) – The mode object containing the electric and magnetic field data to be used for the calculation.- Returns:
The power normalization factor. It will be very close to 1, but it is not exactly 1 due to the numerical integration. currents (dict): A dictionary containing the calculated currents for each polygon. impedances (dict): A dictionary containing the calculated impedances for each polygon.
- Return type:
p0 (complex)
- get_epsilon_rf(frequency, use_charge_transport_data=False, voltage_idx=0)
This function will return the \(\epsilon_{RF}\) tensor with the signature
self.epsilon_rf[vertice_idx, voltage_idx]wherevertice_idxis the index of the vertice in the mesh and voltage_idx is the index of the voltage in the bias points. In case no bias dependent data is available (i.e. no charge transport data is available) thevoltage_idxmust be 0 as the second axis of the array will have size 1.- Parameters:
frequency (
float) – The frequency at which to compute the permittivity tensor. The frequency must be in GHz.use_charge_transport_data (
bool(default:False)) – Whether to use the charge transport data to compute the permittivity tensor. Doing so will yield a \(\sigma(x,y,V)\). Make sure your mesh is appropriate for it.voltage_idx (
int(default:0)) – The index of the voltage to use for the permittivity tensor.
- get_power_loss_per_polygon(mode, frequency)
This function returns the power loss per polygon as described in Pozar [2] eq. 1.92. The result is given in W/cm:
\[P_l = \int_V \frac{\sigma}{2}|E|^2dv + \frac{\omega}{2}\int_V(\Im\{\epsilon\}|E|^2 + \Im\{\mu\}|H|^2)dv\]Warning
This part is not developed yet.
- Parameters:
mode (
Mode) – The mode object containing the electric and magnetic field data to be used for the calculation.frequency (
float) – The frequency at which to compute the power loss. The frequency must be in GHz.
- Returns:
A dictionary containing the power lost per polygon
- Return type:
power_lost_all (dict)
- make_mesh(default_resolution_min=1e-12, default_resolution_max=100, filename=None, gmsh_algorithm=5, global_quad=False, verbose=False, mesh_scaling_factor=1.0)
Returns a gmsh with the geometric information of the photonic device. The mesh generation is currently handled by FEMWELL, but it will later be deprecated in favour of MESHWELL.
- Parameters:
default_resolution_min (
float(default:1e-12)) – The minimum resolution to use for the mesh generation if no resolution is specified for a given entity.default_resolution_max (
float(default:100)) – The maximum resolution to use for the mesh generation if no resolution is specified for a given entity.filename (
str|None(default:None)) – The filename to save the mesh to. If None, the mesh is not saved.gmsh_algorithm (
int(default:5)) – The algorithm to use for the mesh generation.global_quad (
bool(default:False)) – Whether to use global quadrature for the mesh generation.verbose (
bool(default:False)) – Whether to print verbose output during the mesh generation.mesh_scaling_factor (
float(default:1.0)) – The scaling factor to apply to the mesh.
- plot_eps_rf(frequency=10, voltage_idx=0, use_charge_transport_data=False, log_scale_im=True, log_scale_re=True, cmap='jet')
Plots the real and imaginary parts of the RF permittivity tensor.
- Parameters:
frequency (
float(default:10)) – The frequency at which to compute the permittivity tensor. The frequency must be in GHz.voltage_idx (
int(default:0)) – The index of the voltage to use for the permittivity tensor.use_charge_transport_data (
bool(default:False)) – Whether to use the charge transport data to compute the permittivity tensor. Doing so will yield a \(\sigma(x,y,V)\).log_scale_im (
bool(default:True)) – Whether to plot the imaginary part of the permittivity tensor in logarithmic scale.log_scale_re (
bool(default:True)) – Whether to plot the real part of the permittivity tensor in logarithmic scale.cmap (default:
'jet') – The colormap to use for plotting.
- Returns:
The matplotlib figure containing the plots. ax1: The axis for the imaginary part of the permittivity tensor. ax2: The axis for the real part of the permittivity tensor.
- Return type:
fig
- plot_mesh(plot_polygons=True)
Plots the mesh of the photonic device.
- Parameters:
plot_polygons (
bool(default:True)) – Whether to plot the polygons on the mesh. IfTrueit will call theimodulator.RFSimulator.plot_polygons()method.- Returns:
The matplotlib figure. ax: The matplotlib axis.
- Return type:
fig
- plot_mode(mode, Nx=50, Ny=50, xmin=0, xmax=10, ymin=-10, ymax=10, figsize=(10, 4), wspace=0.5, color_polygons='white', color_integral_lines='red', cmap='jet', color_vectors='black')
Plots the electric and magnetic fields of a mode.
- Parameters:
mode (
Mode) – The mode object containing the electric and magnetic field data to be plotted.Nx (
int(default:50)) – Number of grid points along the x-axis.Ny (
int(default:50)) – Number of grid points along the y-axis.xmin (
float(default:0)) – Minimum x-coordinate for the plot.xmax (
float(default:10)) – Maximum x-coordinate for the plot.ymin (
float(default:-10)) – Minimum y-coordinate for the plot.ymax (
float(default:10)) – Maximum y-coordinate for the plot.figsize (
tuple[float,float] (default:(10, 4))) – Figure size as (width, height).wspace (
float(default:0.5)) – The amount of width reserved for blank space between subplots.color_polygons (
str(default:'white')) – Color of the polygons in the plot.color_integral_lines (
str(default:'red')) – Color of the integral lines in the plot.cmap (
str(default:'jet')) – Colormap used for plotting the magnitude of the fields.color_vectors (
str(default:'black')) – Color of the vector field streamlines.
- Returns:
This method does not return any value. It generates a plot.
- Return type:
None
Note
This method uses the mode object to interpolate the electric and magnetic fields onto a rectangular grid.
The fields are split into transverse and longitudinal components for plotting.
The resulting plot includes both the magnitude of the fields and their respective vector field streamlines.
Polygons and integral lines can be overlaid on the plot for additional context.
The projection into a rectangular grid is quite heavy, so it is reccomended to use this method only for small regions to avoid large rectangular grids.
- plot_polygons(color_polygon='black', color_line='green', color_junctions='blue', fig=None, ax=None)
Plots the polygon and line entities in the simulation window.
- Parameters:
color_polygon (default:
'black') – The color of the polygons in the plot.color_line (default:
'green') – The color of the line entities in the plot.color_junctions (default:
'blue') – The color of the junction entities in the plot.fig (default:
None) – The matplotlib figure to plot on. If None, a new figure is created.ax (default:
None) – The matplotlib axis to plot on. If None, a new axis is created.
- refine_mesh(N_nearest_neighbours=200, mode_for_refinement=None)
Refines the mesh based on the computed RF modes.
- Parameters:
N_nearest_neighbours (
int(default:200)) – Number of nearest neighbors when finding the facets in the new mesh corresponding to the old mesh boundaries. If you have a very fine mesh, it is wise to make this number higher. If it’s not high enough, it may not update the named boundaries in its entirety, and you end up with wrong line integrals. Of course, the higher it is, the slower the algorithm.mode_for_refinement (
Mode|None(default:None)) – the mode to be used for refinement.
- Returns:
None
References
R.B. Marks and D.F. Williams. A general waveguide circuit theory. Journal of Research of the National Institute of Standards and Technology, 97(5):533, September 1992. URL: https://nvlpubs.nist.gov/nistpubs/jres/097/jresv97n5p533_A1b.pdf (visited on 2024-05-02), doi:10.6028/jres.097.024.
David M. Pozar. Microwave engineering. Wiley, Hoboken, NJ, 4th ed edition, 2012. ISBN 978-0-470-63155-3. OCLC: ocn714728044.