ECE 5340  FDTD Project:    FINITE-DIFFERENCE TIME-DOMAIN

 

Objectives:

 

¨      Understand and program the FDTD equations in 1D

¨      Observe CW and Pulsed time domain data

¨      Observe numerical dispersion

¨      Understand and program the Mur 1st order absorbing boundary conditions

¨      Understand the relationship between time domain and frequency domain data and use this to calculate reflection coefficient

 

1.     Differential Equations (3D to 1D)

 

                Starting with Maxwell’s equations in the time domain (Ampere’s and Faraday’s laws), differential form, write the 6 coupled differential equations.  (Take the cross products and equate vector components.) 

 

                Convert these equations to the 1-dimensional TE-to-z case by setting d/dy = d/dz = 0 and Ez = 0.  This represents a plane wave propagating in the x-direction.  You should end up with equations for Ey and Hz.  (The TM-to-z case would have similar equations for Ez and Hy.)

 

2.     FDTD Equations (1D TE-to-z case)

 

                Convert the 1D TE differential equations above to their FDTD difference form.  (Use the central difference formula to approximate the derivatives, and solve for Ey ^ (n+1) and Hz ^ (n+1/2).  )

 

                Use the 1D FDTD lattice shown below:

 

 

                                                  |ß----cell I----à|

                                                  |ß----dx-------à|

-------+--------o---------+--------o---------+---------o----------+---------

      Ey(I-1)   Hz(I-1)  Ey(I)    Hz(I)     Ey(I+1)    Hz(I+1)   Ey(I+2) -àI

 

Let the E fields be defined at times n, n+1, n-1, etc.

Let the H fields be defined at times n-1/2, n+1/2, etc.

 

3.     Program the FDTD Equations (1D TE case)

 

                Program the equations in (2) for the geometry shown below.   Use a forced CW source on Ey:  Ey(I=inc) = sin(wt).  Note: although all the test cases here are in air, write your code for arbitrary materials.

 

|--------------------------------|--------------------------------------|

 0                                      Ey(inc)                                                          nx

 

4.     Test the FDTD Equations and observe CW Time Domain Data:

 

                Use F = 2GHz, dx = wavelength/20, dt = dx/(2c), nx=120, inc=60.

 

                Plot the Ey and Hz fields at points A,B,C,D as a function of time for 100 time steps.  (Note that this stops the simulation just before the wave touches the lack of boundary conditions at the ends.)  Give one plot of the four E fields, and another of the four H fields.  Store the E fields at point C for use in problem 7.

 

                Using matlab, make a movie of the E fields at all points as a function of time.

 

                Plot the Ey field at point D against the analytical value:  Ey(x) = sin(wt-bx), where x is the distance from the source.

 

                                                    A   B   C   D

                                                Ey(inc)

|-----------------------------------|----------------------------------|

0                                              inc                                                nx

                                                60                                                120

 

A is located at I=60, at source

B is located at I=63, 3 cells from source

C is located at I=67, 7 cells from source

D is located at I=90, 30 cells from source

 

5.     Observe Pulsed Time Domain Data

 

                Change the source to a raised cosine pulse:

 

                Ey(inc) = 1-cos(wt)  o<t<1/Fmax

                                    0               t > 1/Fmax

 

                Use Fmax = 2GHz, dx = wavelength/20, dt=dx/(2c), nx=120, inc=60.

 

                Plot the Ey fields at points A,B,C,D as a function of time for 100 time steps.  (Notes: If you run more than 120 time steps you will see the waves reflect off the ends of the FDTD mesh.)

 

                Using Matlab, make a movie of the Ey fields as a function of time along the mesh.

 

6.     Observe Numerical Dispersion

 

Use the raised cosine pulsed source, Fmax = 2GHz, dt=dx/(2c), nx=220, inc=110.  Run for 200 time steps.

 

Plot the Ey fields as a function of time 30 cells from the source for dx = wavelength/60, wavelength/20, wavelenth/10, and wavelength/5.

 

Plot the Ey fields 30 cells from the source for the CW source using dx=wavelength/5 and compared to the values observed at point D in part 4.

 

7.     Mur 1st order boundary conditions

 

a)     Make sure you understand the derivation from class:  What approximation is made to make the boundary conditions “first order”?

 

b)     Write the difference form of the 1st order boundary conditions for Ey.

 

c)     Program the 1st order boundary conditions for Ey on both boundaries.

 

d)     Test the boundary conditions:

 

                                Incident Fields:  Use the Ey fields you stored 7 cells from the source in problem 4 as the “incident” fields.  There is no reflection in these fields, because the waves have not yet hit the boundary. 

 

                                Total Fields: Rerun your simulation using the same parameters as problem 4 except nx=30, inc=15, for 100 time steps.  Store the Ey fields 7 cells from the source.  These are the “total” fields, because they contain the incident field plus any fields which might be reflected from the boundaries.

 

                                Reflected Fields:  =   Total fields - Incident Fields

                                Plot the reflected field as a function of time 7 cells from the source.  They should be tiny!  They should look like noise, possibly increasing slightly with time.

 

e)     Run the simulation with boundary conditions for 10,000 time steps to be sure your simulation does not blow up.  Plot the last 3 cycles of the Ey fields 7 cells from the source.  Does it still look like a clean sine wave?

 

8.     Added Source

 

                Program the “added” source in your simulation and test to be sure that it works the same way as the “forced” source in air.

 

9. Numerical Stability

                Test your CW simulations with  several values of dt = dx / K to verify the stability criterion.  For 1D you expect your simulations to become unstable when dt > dx / c

 

10.  Standing Wave

 

                Replace the Mur boundary condition at nx with a metal boundary :

 

                |-------------|------------------------------------||

          Mur        Ey(5) = added source                  Ey(nx=120) = 0.

                                                                            Metal boundary

 

a)     Observe reflected fields:

 

                Using a pulsed source with Fmax = 2GHz, dt = dx/(2c), dx = wavelength/20, plot the Ey fields as a function of time at I=40 and I=100.  You should see the incident and reflected pulses clearly.

 

b)   Verify Propagation Time:

 

                How many time steps will it take for your wave to propagate to the conductor and back to the source using dt=dx/(2c)?_____________  using dt = dx/3c?______________

 

                Verify this using a pulsed source with Fmax = 2GHz, dx = wavelength/20, dt = dx/(2c) and

dt = dx/3c

c)     Observe and understand a standing wave:

 

                Using a CW source with F = 2GHz, dt = dx/(2c), dx = wavelength/20, plot the Ey fields as a function of time at I=40 and I=100.  Indicate the magnitudes (peak value) of these fields on your graph. 

 

                Now have your program find the magnitudes (peak values) of the fields at every point.  Verify the values for I=40 and I=100, compared to those you observed by hand above.  Plot the magnitudes of Ey as  a function of I  (at every point).  You should see the classic ‘standing wave” pattern. 

 

                Where is the first zero of the standing wave in front of the wall? ______

                Where is the first peak? ___________

 

                Plot Ey as a function of time at these points.   This is the total field

 

                Using the method in part 7, find and plot the incident and reflected fields at these two points as a function of time.  Verify that the sum of the incident and reflected fields gives the total fields observed above.

 

11)  Use your program to compute the reflection coefficient from a quarter-wave dielectric transformer as a function of frequency.  A quarter-wave dielectric transformer is a layer of material that can be used to "match" two different materials (such as air and fluid).  The transformer has a characteristic impedance of ht = sqrt(h1h2).  The transformer is a quarter of a wavelength (in the transformer material) long.   Design a transformer to match air (er = 1.0) and water (er = 40.0) at 1 MHz.  The reflection coefficient should be zero at 1 MHz.  Evaluate the reflection coefficient from 0.5 to 2 MHz. 

 

Possible End of Semester Project

 

Research Project:

 

We have talked about the two ways of modeling dielectric interfaces.  Either you can put the boundary of the cell in between the two materials,  and have one cell with each material OR you can average the two materials and put that in for the boundary cell.  The location where you place your source is also critical for this part of the assignment.

 

a)       Draw a very careful cell model of an electric field (transparent source) in front of a dielectric discontinuity, for both possible models.

b)       Write a matlab code to calculate exactly the standing wave magnitude in front of and behind the dielectric interface as a function of distance.   Refer to any basic EM text that includes standing waves.

c)       Compute the magnitude of the standing wave in front of the dielectric discontinuity, and compare it to the analytical solution from (b). 

d)       Determine which of the two models gives you accurate results.  (Hint:  I think that both models give accurate results, if you measure the distance from the interface correctly, and each is different.)

e)       WRITE your results in the form of a short professional report (2-3 pages).  You do not need to explain the FDTD method in detail, but do explain your model in detail, the analytical solution you used, the comparison, and your conclusion.

 

 

 

TURN IN:

 

[ ] Hard copy of your code

 

[ ] Derivations of all formulas used in your code (FDTD, boundary conditions, etc.)

 

[ ] Plots

 

[ ] Answers

 

[ ] Summarize and comment on your results