NMM: Final Project - Modeling of Discrete Electronics

cba

Background

My research involves assembling micro-electronics from pressfit parts. New modeling techniques are required to design and characterize these structures.

The Question

Can a simple and efficient simulation/modeling software be written for electronic digital materials?

My research involves assembling micro-electronics from pressfit parts. Because the structures are inherently three-dimensional, they're not well suited to being designed and modeled in the traditional ways (which, much of the time, assumes relatively planar assemblies. Therefore, new modeling techniques are required to design and characterize these structures. What should these look like?

Who's Done What

The Finite Difference Time Domain (FDTD) method for solving elecromagnetic fields is well studied and quite well understood. The most commonly used FDTD algorithm dates back to 1966 when Kane Yee came up with an efficient and more-accurate formulation for Maxwell's equations in a 3D lattice: Yee's Paper (1966). FDTD is now largely used to study and simulate things like sensitive RF circuits, scattering problems, and metamaterials.

FDTD Solvers:

There are a number of both commercial and open-source FDTD solvers:

References

The two references which I found most helpful (and used significantly in this work are):

The Math - Maxwell's Equations

We'll assume non-magnetic materials ($\mathbf{H} = (1/\mu_0)\mathbf{B}$) and use this formulation of Maxwell's Equations: $$\frac{\partial \mathbf{D}}{\partial t} = \nabla \times \mathbf{H}$$ $$\mathbf{D} = \epsilon_0 \cdot \epsilon_r \cdot \mathbf{E}$$ $$\frac{\partial \mathbf{H}}{\partial t} = -\frac{1}{\mu_0}\nabla \times \mathbf{E}$$ But the electric and magnetic constants differ by many orders of magnitude. This could be problematic when we go to implement it on a computer. $$\epsilon_0 = 8.85 \times 10^8 \:F/m$$ $$\mu_0 = 4\pi \times 10^{-7} \:H/m$$ And so we need a change of variables so that E and H on the same order. Let's try changing E to: $$\tilde{\mathbf{E}} = \sqrt{\frac{\epsilon_0}{\mu_0}} \cdot \mathbf{E}$$ This implies.. $$\tilde{\mathbf{D}} = \sqrt{\frac{1}{\epsilon_0 \cdot \mu_0}} \cdot \mathbf{D}$$ So our update eqations are now: $$\frac{\partial \tilde{\mathbf{D}}}{\partial t} = \frac{1}{\sqrt{\epsilon_0 \mu_0}} \nabla \times \mathbf{H}$$ $$\tilde{\mathbf{D}} = \epsilon_r \cdot \tilde{\mathbf{E}}$$ $$\frac{\partial \mathbf{H}}{\partial t} = -\frac{1}{\sqrt{\epsilon_0 \mu_0}}\nabla \times \tilde{\mathbf{E}}$$ The Courant condition says we don't want electromagentic waves (which travel at the speed of light, $c_0$) to move faster than one unit diagonal per unit time: $$\frac{\Delta{t}}{\Delta{x}} = \frac{1}{\sqrt{n} \cdot c_0}$$ For simplicity I use: $$\frac{\Delta{t}}{\Delta{x}} = \frac{1}{2 \cdot c_0}$$ This makes field updates very simple since $1/ \sqrt{\epsilon_0 \mu_0} = c_0$: $$\sqrt{\frac{1}{\epsilon_0 \cdot \mu_0}} \frac{\Delta{t}}{\Delta{x}} = c_0 \cdot \frac{1}{2\cdot c_0} = \frac{1}{2}$$ So the simple version of the update equation becomes:

curl_h = (hz[i][j][k] - hz[i][j-1][k] - hy[i][j][k] + hy[i][j][k-1]);
dx[i][j][k] = dx[i][j][k] + 0.5*curl_h;
ex[i][j][k] = gax*dx[i][j][k];

Implementation

1D FDTD

I started by implementing one-dimensional FDTD to get a feel for the structure of the algorithm.

1D FDTD, Gaussian point source, no boundary conditions (left), absorbing boundary (right)

2D FDTD

Working my way up, I then implemented the textbook two-dimensional FDTD solver. Here's a plane wave propogating in free space:

For speed I wrote this in C. A python script is passed the name of the C code as an argument and executes bash commands to compile the code, execute it, import the datafile that was generated, and plot or animate it.

3D FDTD

Finally, I made a 3D FDTD solver. While the code takes the same form as the two previous implementations, it's significantly more code and uses a pretty good chunk of memory to run. Once I got the bare bones 3D version working I made a number of additions to make the solver practical for EM simulations.

Perfectly Matched Layer

First, I added a Perfectly Matched Layer (PML) absorbing boundary condition. This treats the outermost cells as lossy material that minimizes the reflection of the oncoming wave.

Total/Scattering Field

Next, I divided the space into two spaces: the total field and the scattered field. The total/scattered field formulation allows a source generated on one end of the domain to be cancelled out on the other end. This minimizes the load on the PML (to minimize reflections) and allows for simulations of things like scattering parameters from plane wave sources.

Here I'm showing a plane wave hitting a dipole antenna:

Graphics - Mayavi

As you can see, I also added the ability to show and animate 3D volumetric data. For this, I'm using Mayavi (and specifically mlab) which can be called straight from the Python wrapper.

Geometry Import - GIF

To make it easy to import geometry into the solver, I added the ability to import volumetric data in the form of animated gifs (or stacks of images). So the GIF on the left gets imported and becomes the volumetric data on the right.

Source:

Disclaimer: I make no guarantees about the working state of these:

Evaluation

To verify the accuracy of the simulation you can either compare the solution to analytical results or to real measurements.

Dielectric Sphere Absorption

I decided to follow the advice of Sullivan (yet again) and shoot a plane wave at a dielctric sphere and compare the EM energy absorbed by the sphere with analytical results (derived from bessel funciton expansions).

I found very good agreement between his results (left) and mine (right)...

S-Parameter Calculation

The next measure I tried, was to calculate the transmission/reflection coeffecients of a pulse propogating down a patch antenna.

Unfortunately, I think there is a bug lurking in my code which is causing a significant amount of ripple on my gaussian pulse. This means that there is a lot of added noise and it is hard to evaluate the accuracy of the reflection.

Implications for Future work

Run more simulations... because it's fun.

Spiral chip antenna:

CBA 2014