Learn how to access and analyze simulation data.
Working with the graphical user interface (GUI):
Understanding analysis and layout modes
FDTD MODE DGTD CHARGE HEAT FEEM INTERCONNECT
This page describes the concept of the Analysis mode. Basically, the software switches to Analysis mode after the simulation has finished running. Typically, analysis mode allows you to view simulation results but not make change to the simulation setup. To change the setup, switch back to Setup mode.
In Analysis mode, simulation object properties can be viewed (
symbol means the object has data or results, similar to havedata and haveresult command) , but not edited. This ensures that at any given time, the simulation data results corresponds to the configuration of the simulation project. Once in the analysis tool, the user continues to analyze the simulation data until they either wish to close the application, or they decide to alter the simulation objects and re-perform the simulation. By exiting the analysis routines and returning to the layout editor, existing simulation data will be erased.
If you would like to keep the existing data and re-run a simulation, the easiest way is to save the simulation with a new file name. Then the new data will be saved to the new file. Another way is to send the data to your Script Workspace, data in this field will not be replaced by re-running the simulation. To do this, right-click the results/data in the Results View , and choose "Send to script".
Simulation files in analysis mode (with simulation results) are typically very large. Switching them back to layout mode makes the file much smaller and easier to share.
Browsing data with the Results Manager
The Results Manager is a tool for viewing simulation data. The Results View window shows all the results for the simulation object that is currently selected in the Object Tree. The Script Workspace works in conjunction with the scripting environment to provide additional GUI-based functionality.
When used in conjunction with the visualizer, the Results Manager provides a very useful and intuitive way of analyzing and visualizing variables and results through the GUI, greatly reducing the need for scripting. The Result View window shows all the results for the simulation object that is currently selected in the Object Tree.
Any simulation object that have results will be displayed with a
symbol on the bottom-right corner (similar to havedata and haveresult command). The name of the available results, and the corresponding dimensions or are displayed in the Dimensions/Value column. One can right click on any of the results to display them in the visualizer, or to send the to the Script Workspace for further post-processing.
With the use of datasets, allowing one to package raw data into meaningful results that can be easily parametrized and visualized. The results for all the standard monitors can be retrieved in the original raw, un-parametrized matrix form (using getdata), or in dataset form (using getresult). For example, in the Result View figure above, the results listed under “rawdata” can be obtained using the “getdata” command. The results listed under "results" are datasets, and can be obtained using the “getresult” command (these calculations will only be carried out when they are visualized). The icons associated with each result reflect the type of result:
String
Matrix : this is a simple matrix result, with no associated parameters
Matrix dataset : this is a parameterized matrix results that contains at least an attribute (result), and an associated parameter
Rectilinear dataset : this is a parameterized matrix result that is associated with a rectilinear grid
Unstructured data : this is a set of data that is not structured in form of a dataset or a matrix and rather consists of several different types of results
For more detail on how to work with datasets in the scripting environment, please see Accessing simulation data. Analysis group objects from the Object library have been updated to return datasets. For an example of how to define dataset results in an analysis group, please see Analysis groups.
The raw data results are all un-parametrized, simple matrix results. To create parametrized matrix datasets from matrices, use the "Send to script" option to copy the variable into the Script Workspace.
Script Workspace
The Script Workspace shows all the variables in the current scripting environment. The variables' current values as well as the corresponding dimensions are shown in a list format. Users can use the visualizer to visualize any variable listed in the Script Workspace by right-clicking on the variable and selecting "Visualize".
The icons in the script workspace above shows that while "sigmaabs" and "sigmascat" are parametrized matrix datasets (since they were the results returned directly by the cross section analysis groups), "sigmaext" is a result that we defined in the script, and is therefore a simple un-parametrized matrix. For example, simply right-clicking on sigmaext and selecting "Visualize" will generate the following plot in the Visualizer:
Here, the extinction cross section is plotted as a function of the index value.
Using the data visualizer and figure windows
FDTD MODE DGTD CHARGE HEAT FEEM INTERCONNECT
This page explains how to use the data visualizer and figure windows to view your simulation data.
Introduction
Simulation results can be visualized using 1D lines, a 2D surface, and 3D plots. These plots can be created from within the results visualizer ("plot in new window" button), or from a script.
Plot windows support standard operations such as axis labels, zoom, export to JPG, etc. Pressing the left mouse button to zoom in by a factor of two, the right button to zoom out by a factor of two, pressing and holding the left-hand mouse button results in a zoom window, and double-clicking either mouse button scales the plot to show all of the data. Some controls are slightly different for each type of plot (see more details in the next sections). For example, only the 3D vector plot provides controls for 3D rotation of the view; SHIFT key + mouse LEFT click allows pan view without 3D rotation. For image plots, the colormap scale can be changed from color to grey, or red2blue.
The FILE menu contains an option to export the current figure to a JPG image file.
The SETTINGS menu contains options for setting the axes, colorbar limits, color map, labels, etc.
Visualizing rectilinear data
In this section, we describe the visualizer tools for data structured in a rectilinear grid (eg. FDTD simulation results).
Available plot types
Line (xy plot): Choose this option to plot a 1D vector versus another 1D vector. For matrices with more than 1 dimension, a slice of the matrix will be automatically chosen. You can use the parameters table on the bottom to choose which dimension to slice and which to plot on the x- or y-axis, e.g., a spectrum (field amplitude vs. wavelength)
Line (Smith chart): Use this option for plotting impedance data. For more information see Wikipedia. (Example: Microstrip with a lumped RLC element)
Line (Polar plot): Choose this option to plot the angular distribution of some quantity. The angular coordinate is in degrees. The data must be provided as a function of the angle in radians. (Example: Mie scattering 3D)
Image: Choose this option to create a 2D image plot of monitor data. For matrices with more than 2 dimensions, slices of the matrix will automatically be chosen. You can use the parameters table on the bottom to choose which dimensions to slice and which to plot on the x or y-axis, e.g., a spatial field profile (field vs. x,y), slice wavelength.
Vector: Choose this option to create a 2D/3D vector plot of monitor data. Choose the parameters table to slice a frequency point. It is possible to create a vector plot of quantities such as the electric field, magnetic field, Poynting vector, etc.
| Note:
|
Features
MATLAB export / MATLAB export script: Only available for line-xy and surface plots. These buttons generate a MATLAB plot using the data plotted in the Visualizer; In addition, the MATLAB export script button creates a script that generates the plot in MATLAB. It requires MATLAB integration. To switch between these two buttons click and hold the button until you see the two options as shown below:
Export to JPG: It exports the figure to a JPG file. For line plots it is possible to export to PDF and EPS files as well. Click and hold the button to see the available options as shown below:
Copy to clipboard: It copies the image to clipboard.
Show/hide chart settings: Opens the Visualizer settings described in the next section.
Get help with graph control options: Available for vector plots only. Shows a description of the control options to zoom in/out, rotate, pan, spin, switch between wireframe and solid view, and reset the camera.
Smith chart additional options:
- ASPECT RATIO: Choose between 1:1 and fill the scene.
- NORMALIZING IMPEDANCE (OHMS): Set the characteristic impedance of the system used for normalization.
Vector plots additional options: Choose to downsample the data, scale vectors by a given factor, show data points, invert background color, apply illumination, and use parallel projection.
Plot data redraw options: There are two options:
- Hold plot settings: When enabled, the current settings (e.g. color bar limits) are applied to any data set (or a slice of it) plotted in the visualizer.
- Auto redraw: When enabled, the plot is refreshed automatically after any change in the settings or data set selection. If the option is disabled, it is necessary to click on the Redraw button to update the plot; this is particularly useful when visualizing large data sets.
Phase shift animation: Visualize a continuous wave response of the field by changing the phase of the propagating field. Works only for scaler data or one of the components (X, Y, Z) of vector data.
- frames per cycle: set number of frames per cycle.
- FPS: set number of frames per second.
- Export: export field with phase changes as a movie.
Export to: Export figure to JPEG, text file or clipboard
Overplot Lines
This option is used to plot a series of lines with different indices from a 2D image plot in a 1D line plot. This is useful to give a clear comparison of the results with respect to different parameter values. To overplot lines:
- Change the "Plot Type" from "Image" to "Line"
- Select the index parameter "Action" to be "Slice"
- Turn on the "Overplot Lines" radio button
- Key in the indices for the lines, separated with a comma
Visualizing unstructured data
The visualizer for unstructured data is very similar to the one for rectilinear data described in the previous section. However, for unstructured data the only plot type available is surface and some toolbar buttons are not available. In addition, the visualizer settings are significantly different as described next.
Some of the visualizer settings are common to 2D and 3D simulation data. These are:
- SHOW: Can be "surface", which will contain the values at each mesh point only, or "surface and mesh", which will superimpose in black the mesh grids on top of the plots, or "mesh only" which will only plot the mesh grid in color.
- AXIS SCALE OPTIONS: Can be "square" which will make sure the two axis are plotted to the same size, or "equal" which will use the same scale for both such that the plot will be to scale.
- LOG SCALE: Will plot the result on a log scale.
- TITLE: The title of the plot can be specified.
- X/Y/Z LABEL: The x, y, and z labels for the plot can be specified.
- COLOR BAR: The color bar limits can be specified and locked. Additionally, it is possible to use a background and border for the color bar.
The specific options for 1D, 2D, and 3D data are described below.
Options for 1D data
1D data captured by line monitors in solvers based on the finite element method can follow a curved path. To express the recorded data in a comprehensive way, the position on the curve is shown as value s normalized to the total length of the curve. Hence, the resulting plot shows the variable of our interest as a function of parameter s that has values between 0 and 1 representing the start and the end of the monitor curve.
The position of the start and end points s=0 and s=1 for closed curves are arbitrarily determined as a unique point with minimal in x, then minimal in y, then minimal in z.
| Note: The exact x,y,z coordinates corresponding to each E field/s data point are available in the monitor results. |
Options for 2D data
It is possible to choose between two plot types: surface and image. In the 2D image plot, the values are indicated using the color bar. In the surface plot, the values are also plotted in the third dimension. For a surface plot, it is possible to plot data as a radiation plot which will plot data in spherical coordinates and the distance of each data point from the origin will determine it's value or intensity. This is specifically useful in the case of far-field projections or a similar type of data where the intensity of data points at each direction in a 3D space is of interest.
Surface (xyz plot) and Image plot of the same data
Image plot and Radiation plot of the same Gaussian beam
Options for 3D data
- DATA VISUALIZATION: Choose between three options: unclipped (default), clipped volume, and clipped plane. When a clipping option is selected a clip plane cuts through the 3D structure and gives insight into the inside of the 3D plot.
- CLIP PLANE: The clip plane can be shown on top of a 3D surface plot. This plane is defined by the coordinates of the origin and components of a vector perpendicular to the plane. For convenience, there are also the quick options for choosing an x-normal / y-normal / z-normal plane. The inside out option will flip on which side of the plane is shown. To change the plane position and/or orientation, move the arrow perpendicular to the plane, grabbing it with a mouse left-click much like a drag action. For Pan view, hold the keyboard SHIFT key and mouse left click.
Clipped volume and Clipped plane
- 2D PLANE DATA: Use this option to export the cross-section 2D data at the clip plane to the script workspace or to show it in a new Visualizer window.
Understanding the Attribute editor
When using line plots, each attribute will appear as a separate line. When using image and vector plots, the selected attribute will be shown.
DATA SET: full data set name (can contain multiple attributes)
- ATTRIBUTE: attribute name
- VECTOR OPERATION: selects a particular component of a vector attribute e.g. (Ex, Ey, |E|^2)
- SCALAR OPERATION: selects a particular component of a scalar attribute e.g. (real, imag, abs, angle)
- SCALE: multiplier for the data being plotted
- LEGEND: this name will be shown in the legend of the plot
- NOTES: additional information added by the user about the attribute
- VIEW DATA: allows users to view the data in a table format as shown below
In this table format, users can select any portion of the data and "Copy" or "Export" it into a text file. Alternatively, users can also send any portion of the data into the Script Workspace.
Understanding the Parameter window
- ATTRIBUTES: Name of the associated attribute
- PARAMETERS: Name of the parameter
- ACTION: Control how the parameter is treated in the plot. For example, select which axis to plot the parameter on.
- VALUE: displays the value if it is a singular value or is blank if there is a vector of values
Plotting multi-dimensional attributes
Line plots
One parameter must be selected to plot on the x-axis. All other non-singleton parameters must be set to Slice. A specific slice can then be selected for those parameters.
Surface plots
Two parameters must be selected to plot on the x- and y-axis. All other non-singleton parameters must be set to Slice. A specific slice can then be selected for those parameters.
Vector plots
The spatial dimensions (x,y,z) are always selected to plot on the x-, y-, z-axes. All non-spatial parameters must be set to Slice. A specific slice can then be selected for those parameters.
Analysis groups
Calculating charge distributions and currents in metals induced by an optical field
The following pages describe how to make charge distributions and current measurements within metals. Two simple examples are provided below.
| Note: This analysis should only be applied to real metals (i.e. materials where it is reasonable to assume the material is a cloud of free electrons). This analysis will not work for PEC materials. |
Background
To calculate the charge distributions and current densities, we treat each metal as a cloud of free electrons, i.e. a plasma. To calculate the current density in a plasma we first recognize that all material properties within the FDTD simulation are implemented via an effective material permittivity:
D=εmaterialE
For the purposes of this calculation, we assume there are two contributions to the material permittivity: the background permittivity (in this case, the permittivity of free space) plus the contribution from the current density.
D=(ε0+εplasma)E
The first term is Do, the displacement field that would exist in free space for the given electric field. The second term is proportional to the current J. For more details, see the notes below.
D=εmaterialE=ε0E+εplasmaE=D0+iJω
Solving the above equation for J, we get
J=−iω(D−D0)=−iω(εmaterial−ε0)E
Note: The total material permittivity is created from two contributions: one from the polarization of the medium due to bound charge, and one from the current density due to free charge. Ampere's law can be written as
−iωDbackground=∇×H−J
Dbackground=ε0E+P=ε0(1+χ)E
We can rewrite this equation as
−iωDbackground=∇×H−J−iω(Dbackground+iJω)=∇×H−iωDmaterial=∇×H
and, assuming that J is proportional to E, we have
Dmaterial=ε0E+P+iJω=(ε0+ε0χ+εplasma)E=εmaterialE
For plasma materials, such as metals, we have assumed that the susceptibility of the medium is 0 (i.e. P=0).
Charge and current in a nanoparticle
This page describes how to calculate the charge and current density inside metals. The associated simulation files use the current charge density analysis object for calculating these quantities. This analysis object can be found in the object library in the Advanced analysis section.
Example
The associated files section provide a 2D and 3D example simulation using the current charge density analysis groups. You can copy and paste these groups into your simulations or insert them directly from the object library.
| Note : Mesh size and mesh refinement setting It is often necessary to use a small mesh size for simulation involving metals, particularly when the metals have rounded edges. The Conformal meshing algorithm feature can provide improved accuracy at interfaces, but it can also introduce additional errors into the simulation when applied to metal interfaces. Therefore, the default settings in FDTD is to apply the Conformal mesh to dielectric interfaces, but not to metal interfaces. In these particular example files, the mesh refinement option has been set to 'Conformal Variant 1' (CV1), rather than the default setting of 'Conformal Variant 0' (CV0), meaning that the conformal meshing algorithm will be applied at the metal-dielectric interfaces. This change should only be made to your simulation file after careful convergence testing to confirm that the CV1 option gives better (i.e. more converged) results than CV0. Until you do this type of convergence testing, we recommend using the default CV0 option in your simulation. The factors that make CV1 appropriate for this simulation include:
For more information on the mesh refinement options, see Mesh refinement options . |
To reproduce the following figures, open and run usr_current_density_2D.fsp . After the simulation is complete, you can reproduce the following figures by visualizing the object data. The following results show the current and electric field data for lambda = 340nm. Also note that the image colorbar scales have been adjusted.
|
|
| The X and Y components of the current density J. These figures show the current flowing back and forth in the X direction as the wave propagates along the particle in the Y direction. Also notice that the current is zero outside the particle, as expected. |
|
|
| The X and Y components of the Electric field. The current density is proportional to the electric field, except for a 90 degree phase shift. |
The same analysis group also calculates the divergence of the electric field, which gives the total charge density. This is shown in the plot below, we can see that the charges are confined to the edge of the metal, as expected. The free charge density is calculated using rho = eps0*|div(E)|. This assumes that the material is metal with no bound charges, hence eps0 is used instead of eps_material.
| Note : Using symmetric/anti-symmetric boundary condition If symmetric/anti-symmetric boundary condition, make sure the symmetry condition (Analysis tab->Script tab) is specified in the analysis group to avoid matrix size mismatch between parameters used. |
Current in a wire
This page describes how to calculate the total current within a metal wire. The associated simulation file uses the induced current analysis group that can be added to your file from the Advanced Analysis section of the object library.
Example
The analysis group in the associated file uses Ampere's law to calculate the current induced in a metal wire by the incoming EM fields.
I=∫LH⋅dL−d(∫SD⋅dS)dt
| Note : This analysis object treats all metals as plasma materials, with a refractive index of 1. For the purposes of this calculation, metals are defined as materials with absorption. This generally works well when the background material is a non-dispersive dielectric material (i.e. air or glass). However, if your background material is a dispersive (lossy) dielectric, the analysis script will need to be modified. |
To run this example, open the attached simulation file and enter the following script commands:
run;
runanalysis;
f=getdata("induced_current","f");
current=getdata("induced_current","I");
plot(f/1e9,abs(current),"f (GHz)","|current|","Induced current");
Creating curved and angled monitors
This page shows two examples to explain the idea of creating curved or angled monitors using groups and interpolation. This approach is required because the basic monitor objects can not be curved or rotated.
Example 1 (using groups)
This example shows how to use groups of point monitors to obtain frequency domain data along arbitrary paths and arbitrarily oriented planes. In the example below, the mode profiles for a bent wire waveguide are determined at 0 degrees and at 45 degrees into a bent waveguide. In addition, the electric field intensity is determined along a circular arc in these planes.
To obtain the images shown below, run the usr_curved_monitor_group.fsp simulation. Then, run the usr_curved_monitor_group.lsf script file. Please note that to save on meshing time, the mesh refinement setting was set to staircase. The first two images shown below come from the two monitor groups located at the input section of the waveguide. The left image is obtained from the monitor plane2 group. Since the input to the waveguide is normal to the X-axis, the same distribution can be obtained with a conventional 2D X-normal power monitor. On the right, you can see the electric field intensity obtained from monitors that form the circular arc which is drawn in red in the left image. These monitors are contained in the curve2 group, and the location of the monitors can be easily modified to lie along any arbitrary curve in 3D.
The next two images shown below come from the monitor groups located at 45 degrees into the waveguide. In this case, the 2D plane monitor data cannot be obtained from a conventional monitor. As with the images above, the image to the right shows the electric field intensity around the red circular arc from the left image.
The image below (left) shows the transmission vs. frequency for both the monitor group (composed of 50x50 profile point monitors) located at 0 degrees and the 2D X-normal power monitor. The transmission from the two types of monitors will converge to the same result as the number of monitors in the monitor group increase. The right image shows the same plot for the monitor group at 45 degrees.
Example 2 (using interpolation)
In some situations, it is useful to know the field profile or other data along an arbitrary path. Unfortunately, 1D line monitors in FDTD/propagator are always straight, and always parallel to the X,Y or Z axes. If you need data along some other path, one option is to use a 2D monitor. When the simulation is finished, the interpolate function can be used to get the data along an arbitrary path within the plane of the 2D monitor. In the following example, we get the field profile along a straight line with a 45 degree angle, and a semi-circle. Interpolation along a curve can be done using the interpolation analysis object from the Data representation section of the object library.
Interpolating monitor data along an arbitrary path.
The data interpolation can either be done using a script file or with the interpolation analysis object. After running the associated simulation file, you can either run analysis using the interpolation object to get the plots, or else run the associated script file. The script first plots the full 2D field profile that was recorded by the monitor (shown above). Next, we define a new set of x,y position vectors that define the path of interest. The first example is a line from (-5,-5um) to (5,5um). The interpolation function is used to get the field values along this line. When the interpolation is finished, the path and the field data along that path are plotted. Finally the calculation is repeated for a curved path along a semi-circle with a center of (1,0) with a radius of 4um. The two sets of interpolated data are shown below.
| Note: Minimizing data collected If you need to minimize the simulation time or the amount of data collected, this method of recording 2D data then interpolating to a 1D line is not ideal. An alternate approach is to create a line of point monitors located directly along the path of interest. See the previous example for details. |
Tip for adding structure outlines to field plots
This page provides an example file which uses the E-field outlined analysis object that creates a plot of |E|^2 with the structure outline. The E-field outlined object can be inserted from the object library in the data representation section. The files in this example were created using FDTD, but the same analysis group can be found in the Component library in MODE.
If you wish to make more sophisticated structure outlines on your figures, the best option is to export the data to another data visualization tool such as Matlab. In addition to exporting the field data, you must export the refractive index data. The refractive index profile can be used to create the structure outline. In Matlab, you might use the contour functions.
The associated simulation file has an analysis group that returns the field profile with the structure outline super-imposed on the data. To reproduce the above figure, open the simulation file and run the associated script. It will create the following figures: |E|^2 with and without the structure outline.
Changing the far field refractive index analysis object
By default, far field projections assume that the material at the monitor location extends to infinity. In the following figure, this implies the substrate material extends to infinity. Obviously this is not always true.
This page describes how to calculate the far field distribution assuming the refractive index in the far field is different from the index in the near field, which makes it possible to include effects such as the Substrate-Air interface shown above.
Note: The descriptions and examples of the far field projection calculation on the following pages are primarily intended for users of FDTD. For users interested in calculating far field projections with MODE, these descriptions are basically still correct, although some subtle differences do exist.
Two methods are available: Directly setting the far field refractive index in the far field projection functions, and applying the Fresnel equations to the far field data in an additional post processing step.
Setting the far field refractive index in the far field projection functions
The far field projection functions have an optional argument to set the far field refractive index:
out = farfield3d("mname",f, na, nb, illumination,
periodsa, periodsb, index, direction);
To understand this option, it is important to remember that the far field projection calculation is basically a plane wave expansion of the near field data. To calculate the expansion, it is necessary to know the background refractive index, since the ratio between E and H fields in a plane wave depends on the refractive index. By default, the refractive index at the monitor location is used for the expansion, but the 'index' property of the projection function allows a different index to be specified.
Expanding the fields using a different refractive index can provide useful information, but some aspects of the data (particularly the field amplitude and power measurements) can be slightly challenging to interpret. To illustrate the basic issue, consider a single plane wave in air, propagating in the forward direction. If this plane wave has an electric field amplitude of 1 V/m, the magnetic field amplitude will be sqrt(eps0/mu0) = 0.0026.
Now imagine expanding this field profile, using plane waves that exist in a medium with a refractive index = 2. The key point is that the ratio between E and H will be different in this medium by a factor of the the refractive index: If E=1, then H=2*sqrt(eps0/mu0) in this medium. To represent the original field profile using these plane waves, it is necessary to combine a forward propagating wave with a electric field amplitude of 0.75 propagating in the forward direction and a backwards wave propagating with an amplitude of 0.25. Combining these two plane waves allows us to reconstruct the original field profile. Notice that the sign of H is reversed for the backward propagating wave.
E field: 0.75 + 0.25 = 1
H field: 0.75*2 - 0.25*2 = 1
This method can be used to quickly see how refraction from a far field interface will affect the angular distribution of radiation. It is also used in various types of advanced data analysis. However, it is generally not the best option for calculating the effect of a far field interface, especially when the results will be compared to experimental measurements (i.e. power or field intensity measurements). Instead, the Fresnel correction method described below should be used.
Example
The associated simulation file has a gaussian beam propagating at a 10 degree angle of incidence in a medium with a refractive index of 2. The analysis script will plot the far field for a refractive index of 2 and 1.
Efar = farfield3d("T",1,res,res);
Efar_air = farfield3d("T",1,res,res,1,1,1,n_air);
Notice how the angle of the beam changes. In the FDTD simulation (with a refractive index of 2), the gaussian beam propagates at an angle of 10 degrees. The standard far field projection in the substrate shows the beam continues to propagate at a 10 degree angle. When the far field refractive index is changed to 1, the angle of the beam shifts to about 20 degrees. Snell's law can be used to confirm this is the expected shift:
θ2=asin(n1sin(θ1)n2)=asin(2sin(10)1)=20.3deg
Fresnel correction
Alternatively, it is possible to calculate the far field data using the 'near' field refractive index, then use Snell's law and the Fresnel equations to account for the far field interface. Snells' law is used to calculate the change in propagation direction that occurs when the fields pass from one material to the other. The Fresnel equations are used to account for the fraction of power that is lost due to reflection from the interface. The details of this calculation can be found in the associated example files. In particular, see the analysis script of the 'far_field_change_index' analysis group.
This approach can be slightly more complicated that the first method, but it is recommended when it is necessary to get correctly normalized field amplitudes and power measurements. The example script will calculate the |E|^2 in the far field using the default projection (assuming the substrate with an index of 2 extends to infinity), using the plane wave expansion method with a far field index of 1, and using the Fresnel correction using a far field index of 1.
| Note: It is important to understand that multiple reflections (between the interface and the device in the FDTD simulation) effects are not taken into account by this technique. Fortunately, such reflections are often small, making this approximation valid in many situations. |
Far field projections from a box of monitors
This page describes how to calculate fields outside of a closed surface using a box of monitors and the far field projection functions.
Lumerical provides many built in analysis groups in our object library. Please press this button
to open the online library of analysis groups and select the far field category to see which analysis groups are available.
| Note: The descriptions and examples of the far field projection calculation on the following pages are primarily intended for users of FDTD. For users interested in calculating far field projections with MODE, these descriptions are basically still correct, although some subtle differences do exist. |
Using the surface equivalence theorem, it is possible to show that fields radiated outside of a closed box by sources located inside the box can be determined exactly from the field components at the surface of the box. Since Maxwell's equations are linear, the fields outside of the box can be computed by calculating the far field projections for each surface of the monitor box and then summing the results.
The associated files for this section, the usr_farfield_angular.lsf script file and usr_farfield_symmetry.fsp simulation show how to use this method to obtain the electric field and intensity in the far field due to two point sources. For sake of simplicity this simulation only contains two dipoles surrounded by a box of monitors. However, this method also works when there are structures located in the monitor box. See, for example, the Mie 3D application example.
A perspective view of the simulation is shown below to the left. The yellow box in the image shows the edges of the monitors, which are grouped together in an analysis group. When the simulation has been run, the usr_farfield_angular.lsf script file can be used to run the analysis script and create the plot of the far field data shown below to the right. All of the far field projection script is contained in the analysis group, and the script file is used to run the analysis script plot the resulting data.
| Note: Particle scattering with a substrate The 'box of monitor' far field projection is frequently used when studying particle scattering. This is fine, but it's important to remember the requirement that everything beyond the monitor box must be a single homogeneous material. When a substrate is present, this approach is not valid. Instead, the best way to calculate the far field scattering pattern is to use a single monitor, located above or below the particle (depending if you want scattering in the forward or backwards direction). You can then use the standard farfield3d function. When using a single monitor, it's important to make the simulation span large enough that most of the scattered light will pass through the monitor before hitting the PML absorbing boundary conditions. |
| Note: Negative signs in front of far field projections when the normal to the monitor box points along the negative axis. In the case where no symmetry is used, the scat analysis group contained in the usr_farfield_symmetry.fsp simulation computes the far field projection from the sum of the far field projections from the six monitors which make up the monitor box. Notice the negative sign on some terms. E2_far = farfieldexact("x2",x,y,z) - farfieldexact("x1",x,y,z)
+ farfieldexact("y2",x,y,z) - farfieldexact("y1",x,y,z)
+ farfieldexact("z2",x,y,z) - farfieldexact("z1",x,y,z); The negative sign is required because the projection needs information about the projection surface normal. In the image below, the projection normal points outwards from the monitor box for the left monitor (x1) and right monitor (x2) are depicted with red arrows. However, the far field projection function does not know which monitor is on the left; it always assumes the monitor normal is pointing along the positive axis, as depicted with the yellow arrows. The negative sign in front of the far field projection for the x1 monitor in the equation above corrects for the fact that the yellow and red arrows point in opposite directions.
Only the farfieldexact, farfieldexact2d and farfieldexact3d commands allow projections from multiple monitors added to create a total far field projection. |
| Note: Far field projections and symmetric/anti-symmetric FDTD region boundaries When symmetric or anti-symmetric boundary conditions are used in FDTD simulations, it is possible that a monitor will lie entirely outside of the simulation region. Although the monitor does not collect any data, it is possible to obtain the monitor data from other monitors in the simulation. In the usr_farfield_symmetry.fsp simulation, there are two monitors which lie outside of the simulation region, and the analysis object uses data from other monitors to obtain the contribution from the monitors which are not included in the simulation. |
| Note: Far field half space analysis vs resolution and number of frequency points This calculation will take longer than the polar plot. In particular, the resolution and the number of frequency points can significantly affect the time it takes to computer the half space. Reduce either or both parameters can noticeable speed up the analysis. |
Related publications
Allen Taflove, Computational Electromagnetics: The Finite-Difference Time-Domain Method. Boston: Artech House, (2005).
Far field directivity calculations of an antenna
This page describes how to calculate the directivity of an antenna, defined by IEEE as the ratio of the radiation intensity in a given direction from the antenna to the radiation intensity averaged over all direction, in FDTD by using a set of monitors and the far field projection functions built into the directivity analysis group. The user will learn how to obtain the azimuth and elevation directivity patterns, radiated power, and maximum directivity of an "ideal" dipole source using the directivity analysis group.
Projections onto a Spherical Surface
The directivity analysis group is a natural extension of the Projections from a monitor box analysis group in which surface equivalence is once again invoked to determine the far fields radiated outside of a closed box from the radiating sources (Es , Hs) located inside the box, as depicted in the above figure. These farfields (Er , Eθ , Eϕ) are projected onto a sphere of a sufficiently large radius (in this case, 100meters) using the farfieldexact script command and related to the angular distribution of radiation intensity (U) - the power radiated from an antenna per unit solid angle.This allows us to determine the directivity (D) which in mathematical terms is written as D=4πU/Prad
in which Prad is the power radiated out of the monitor box. Often the maximum directivity is reported for an antenna, which is defined as maximum (Dmax) over all angles.
Infinite Ground Plane
Illustration of a set of monitors above an infinite ground plane and their virtual images.
In many instances, it is appropriate to consider the antenna as mounted on top of a metallic ground plane. When the ground plane is sufficiently large and we are operating at low enough frequencies, it may be modeled as an infinitely large perfect electrical conductor (PEC). The presence of the PEC ground predictably alters the radiation properties of the antenna system. For instance, any energy radiating from the antenna directed towards the ground plane undergoes a reflection and reflects back, and the antenna will setup currents on the ground plane which radiate themselves.
Image theory is used in the directivity analysis group to account for the infinite PEC ground plane. The antenna over the ground plane is replaced with an equivalent system of the antenna and its virtual image, which provides the correct field distribution above the ground plane. Below the ground plane, the fields are zero. This technique is achieved in the analysis group by projecting the set of monitors above the ground plane into the negative half-space (-z), which is equivalent to having a set of monitors below the ground plane and projecting them into the positive half space (+z). This procedure allows us to use a finite-sized box to capture all the radiated fields, even though induced currents on the infinite PEC ground plane contribute to the radiation outside of the box.
User Input and Results
After inserting the directivity analysis group, a set of monitors forming a box will be inserted into the simulation domain. In the Setup Tab of the analysis group, the user has must specify the span of the monitor (x span, y span, z span) and the down sampling of the fields on the monitors. Keeping down sample at 1 results in no down sampling. If an infinite ground plane is in the simulation, the variable inf gp should be set to 1.
| Note: Monitor and Antenna Setup For Infinite Ground Plane When the user specifies an infinite ground plane, the analysis group assumes the ground plane’s surface normal is in the z direction and is located on the z=0 plane. Furthermore, the z span variable specifies the span of the monitor starting from z=0 to z=z span. |
If the power from the antenna’s feeding source passes through the box monitor and is not to be included in the directivity calculation, the window source variable should be set to 1. The user must specify which monitor the window should be located on by setting the window surface variable to x1, x2, y1, y2, z1, and z2. The size and center position of the monitor (window center (i, j), window width, window height) fix the size and center position of the window. For example, if window surface=z1, window center (i,j) should be entered as the x-axis position and y-axis center position of the window and window width and height as the x span and y span of the window. The variables coordinate entries (i,j) and the window's surface normal form a right handed triplet.
Under the Analysis Tab of the directivity analysis group, the user can specify the far-field resolution in θ and ϕ (in degrees) and the frequency to obtain the far-field results at.
| Note: Directivity in Freespace In many antenna applications, the directivity is measured and/or calculated in free space and not inside other dielectric mediums. The directivity analysis group assumes that the radiation is occurring into either unbounded free space or, in the case of the infinite ground plane, a bounded half-space. |
| Note: Resolution vs analysis time The resolution can significantly affect the time it takes to compute the directivity. Reducing the resolution can noticeable speed up the analysis, but can have a significant impact on the calculated radiated power and the absolute value of the directivity. Care should be taken to make sure any fine features in the far field are properly resolved with the specified resolution. |
The results view of the analysis group provides the θ (Dθ) and ϕ (Dϕ) polarizations of the directivity in linear units and the total radiated power (Prad). The total directivity (Dtot) is simply the sum of the two polarizations.
Radiated Power
The calculation of Prad assumes that all the power was captured by the box of monitors and that this captured power is from the radiating antenna element itself and not from external sources from outside the box. However, in practice the antenna will need to be excited using an external source (termed an antenna feed), which must pass through the box. The power carried by the antenna feed will significantly influence the calculated Prad. To correct for this, we allow the user to specify a window around the antenna’s feed which subtracts out the power carried by that feed line. The size of this windows must be carefully chosen to capture only the input power and not the radiated power from the antenna. In most cases, the location and size of the source window and source should match (see the Rectangular Probe Antenna example for additional details). In cases where the antenna feed itself contributes to the total radiation, the impact of the window's size on the directivity and radiated power should be tested.
| Note: Radiated Power (Prad) Normalization The calculation of radiated power is similar to the approach used in the Power transmission box analysis group. The transmission is measured through each monitor (subtracting our the window if needed) and then the radiated power is calculated with the formula Prad=trans×sourcepower |
| . When using the mode source or ports, the input power into a waveguide or transmission line may not equal the sourcepower which can introduce a normalization error into Prad. One way to correct for this error is to re normalize the radiated power to the source's/port's input power, as done in the Quarter-wave Monopole example. |
Ideal Dipole Source
In Projections from a monitor box , the far-field patterns are found for two dipoles surrounded by the box monitors. Here, we use the directivity analysis group to investigate the radiation performance of a single dipole radiating into unbounded free space and a single dipole radiating into free space bounded by a PEC infinite ground plane for z ≤ 0. The associated files for this section, dipole_directivity.lsf and dipole_directivity.lms are used to find the directivity patterns and radiated power of the dipole antenna. A perspective view of the two simulation models is shown below. The yellow box shows the edges of the monitors, which are grouped together in the analysis group. Whereas the image on the left for the unbounded antenna reveals the box completely encloses the antenna, the images on the right for the bounded antenna reveals the box only partially encloses the antenna, in which the bottom edges of the monitors touch the PEC ground plane.
Running the script runs the simulation file for both bounded and unbounded simulation setups and generates the plot results shown below. It is evident there is a strong match between the theoretical (derived from [1] Balanis) and the FDTD results. Although only Dθ is shown, Dϕ from FDTD is well below than -120 dB, which is well within numerical error to the theoretical value of −∞
dB.
=======Unbounded Dipole========= Incident Power: 10.015 fW Radiated Power: 9.996 fW Radiation Efficiency: 99.8 % Directivity: 1.761 dB
=========Bounded Dipole========= Incident Power: 5.008 fW Radiated Power: 4.998 fW Radiation Efficiency: 99.8 % Directivity: 4.771 dB
The script also generates the radiation performance of the unbounded and bounded antenna. As expected, in the absence of loss the radiated and incident power are almost equal, resulting in a radiation efficiency of nearly 100%. Since the dipole source is used, the incident power is calculated from the dipolepower script command at λ=0.45μm. The calculated directivity of the unbounded and bounded antennas is extremely close to their theoretical directivity of 1.761dB and 4.761dB, respectively.
Related publications
[1] C. A. Balanis, Antenna Theory and Design, 4th Edition. John Wiley & Sons (2016).
Calculating the far field polarization ellipse
This page provides an example which uses the polarization ellipse analysis object that creates a plot of the polarization ellipse of a grating order. From the ellipse, it is easy to determine the primary polarization angle and the degree of circular polarization. In addition, this analysis group output the polarization handedness for all grating orders and wavelengths. In this example, the point of view of the receiver (against the propagation direction) polarization convention is followed.
The polarization ellipse analysis object can be found in the object library
in the far field projections section.
Simulation setup
Structures
We create a nanowire polarizer with a series of rectangles of Perfect Electrical Conductor (PEC). These wires reflect S polarization, but allow P polarization to pass through. It is always a good idea to start an investigation with a very simple test case. Before studying the nanowire polarizer, we will check our analysis script by simulating a plane wave propagating through empty space. This is an ideal test case because we know the measured polarization should match the source polarization. Note, the structures outside the simulation region are just for visualization purpose.
Sources
The simulation uses two plane wave sources which allow us to simulate linear, circular, or elliptically polarized light in a single simulation. For example, to create circularly polarized light, set the amplitude of each source to 1. Set a 90 degree phase difference between the two sources. See the Circular polarization page for more information.
Simulation region
We use Bloch boundary conditions because the source is at a non-normal incidence.
Monitors
The Polarization ellipse analysis object is composed of a power monitor and some associated script commands. Lumerical provides many built-in analysis groups in our object library. Please press this button
to open the online library of analysis groups and select the far field category to add a similar object to any simulation.
Results - polarization ellipse
It is always a good idea to start with a simple test case to ensure there are no errors in the simulation setup or analysis scripts. Therefore, we will start by disabling all of the structures in the simulation. With no structures, we can expect the polarization measured by the monitor to match the source polarization. By changing the source amplitude and phase, you should be able to reproduce the following results. Since there is no structure in this simulation, all of the light will be transmitted to the zeroth grating order in far field.
The title of the plot should indicate the grating order of interest, the grating order angles, the ellipse major axis angle and the major/minor axis ratio.
Disabling all of the structures (test case)
|
| Linear polarization, 0 degree (P polarization) Source settings:
|
|
| Linear polarization, 45 degree Source settings:
|
|
| Circular polarization Source settings:
|
The above tests indicate that our simulation and analysis script are correctly setup. We can now move on to a slightly more complicated system, with the nanowire grating structure. This structure acts like a polarizer, allowing the P polarization to pass through, while reflecting most of the S polarization. After running the simulation with a circularly polarized source, we get the following polarization ellipse. The ellipse clearly indicates that most of the S polarization has been removed, as expected.
Enabling the structure
|
| Circular polarization Source settings:
|
Extracting Gs, Gp data to text file
Since the analysis group only plots the polarization ellipse for 1 grating order at a time, one may want to extract all Gs, Gp information as a function of all grating orders and frequency points. User may use the below example command to export Gs_plot, and similarly for Gp_plot, to a text file (and then to Excel file for further analysis). Note that, the write command cannot export multi-fold matrix data to a single .txt file. The extraction of the first (Gs), second (n) and third (m) dimensions data is shown, user may further modify these lines to extract the fourth dimension (frequency). Please see, for , write , pinch , text file .
#clear; # may need to clear the script workspace
Gs_plot=getresult("polarztn_ellipse","Gs_plot");
Gs_plot=pinch(Gs_plot.Gs_all,4); # only pick the first frequency point
sz=size(Gs_plot); # size of each dimension
for(i=1:length(sz)) { # loop over the last dimension
write("Gs_plot.txt","slice "+num2str(i)); # slicing the third dimension
write("Gs_plot.txt",num2str(pinch(Gs_plot,length(sz),i)));
}
Results - Grating order strength and polarization handedness
The Grating order strength figure indicates that a fraction of the power is scattered to the (-1,0) grating order (about 60%). To check the polarization of this order, visualize pol in the Result View, then choose G to be the attribute. To see the exact numbers, press View Data in the attributes section. It is not necessary to re-run the actual simulation. To slice it over other frequency points, choose "lambda" in the parameter section and slice it to your wavelength of interest. In this example, the point of view of the receiver (against the propagation direction) polarization convention is followed. For detailed definition of the output quantities, please refer to the analysis group, Analysis -> Script tab.
Absorption per unit volume
This section describes a number of techniques for calculating power absorption per unit volume (Pabs) due to material absorption. The files in this section were created using FDTD, but a similar approach can be applied for MODE Solution's propagator. Once Pabs is known, it is easy to calculate absorption in some volume by integrating the loss function over that volume.
Lumerical provides many built in analysis groups in our object library. Please press this button
to open the online library of analysis groups and select optical power category to insert the analysis groups shown here.
The following sections provide three techniques for calculating power absorption per unit volume (Pabs).
Simple method
This section describes the standard technique used to calculate power absorption per unit volume. Power absorption is proportional to the electric field intensity and the imaginary part of the permittivity. Therefore, the calculation requires field data from a frequency domain power monitor, and permittivity data from an index monitor.
Advanced method
The advanced method is an extension of the simple calculation that minimizes interpolation errors near interfaces. This technique is appropriate for simulations with large field discontinuities at the structure interfaces. Such discontinuities are common in simulations involving metals.
Divergence of Poynting vector
This section describes how to calculate the absorption per unit volume from the divergence of the Poynting vector. Divergence calculations tend to be very sensitive to numerical problems, which makes this technique inferior to the other techniques described in this section. This method is presented only because it was a recommended approach in earlier versions of FDTD.
Calculating absorbed optical power - Simple method
This page describes a simple and easy to understand approach to calculate absorbed optical power per unit volume (Pabs) due to material absorption. Once Pabs is known, it is easy to calculate total absorption in some volume by integrating the loss function over that volume. A more numerically accurate version of the calculation is described in Advanced method. Calculating the divergence of Poynting vector is an alternate approach to get Pabs.
Theory
The absorption per unit volume can be calculated from the divergence of the Poynting vector,
Pabs=−0.5 real(∇⋅P)
It is possible to calculate the absorption directly from this formula (see the Divergence of Poynting vector section), but divergence calculations tend to be very sensitive to numerical problems. Fortunately, there is a more numerically stable form. It can be shown that the above formula is equivalent to
Pabs=−0.5 real(iωE⋅D)
With a little more work, we get the desired result
Pabs=−0.5 ω|E|2 imag(ε)
To calculate the absorption as a function of space and frequency, we only need to know the electric field intensity and the imaginary part of the permittivity. Both quantities are easy to measure in an FDTD simulation.
Simulation setup
The example simulation is composed of an Air - Glass - Si stack. A focused beam is incident on the stack. Air and Glass are not lossy materials, so we only expect absorption to occur in the Si layer.
A 3D power monitor is used to measure the electric field intensity. Similarly, an index monitor is used to measure the permittivity. Both monitors should be the same size. 3D monitors tend to require large amounts of memory, so they should be kept as small as possible. In this example, both monitors cover the entire simulation region, but this is only to demonstrate that the loss is zero in the Air and Glass regions.
The simulation setup is shown in the following figure.
| Note: Check memory requirements Since 3D monitors require large amounts of memory, be sure to use the Check memory requirements button before running the simulation. In this example, about 100MB of data is collected. Up to 1GB is required to run the analysis script. Memory requirements for larger simulations can be much higher. |
| Note: Reducing memory requirements If the memory requirements are too large, consider making some of the following modifications:
|
Calculating Pabs
Run the simulation with FDTD mesh accuracy of 2 and then run the analysis script, which will create a number of figures. The first three figures show a cross section of the electric field intensity and material properties at y=0 for the second frequency point (at a wavelength of 500nm).
|
| The electric field intensity. The figure can be divided into 4 regions.
|
|
| The real part of the refractive index. There are three materials in this simulation.
|
|
| The imaginary part of the refractive index. The imaginary part is a measure of the material loss. Loss will only occur in materials with k>0.
|
From this data, the script calculates the material absorption as a function of space and frequency. A cross section of the absorption function is shown below. Notice that loss only occurs in the Si layer.
Total absorbed power
The total absorbed power can be calculated in two ways:
- Integrate the loss function over the entire simulation volume.
- Use the transmission function and a set of 2D power monitors to measure the amount of power flowing into and out of the Si layer. The difference between these values will be the material absorption. In this example, T2 measures the power flowing into the Si, and T1 measures the power flowing out. To be strictly correct, there should also be monitors at the sides of the Si layer, but they were not added to this simulation because most of the power is flowing in the Z direction.
Both techniques give similar results. If a smaller mesh size was used (default mesh accuracy=2), the results would be even closer. Also, although we don't see it in this example, the discrepancy can increase slightly at the shorter wavelength where the absorption of the Si is higher. This is due to some numerical errors in integrating the loss in the first dz of the Si and can be improved by using the Advanced method.
Absorbed power in spherical region
When calculating absorption within a rectangular region, a box of 2D power monitors is computationally more efficient (since it only requires surface data rather than volume data, and because the index monitor is not required). The primary reason for calculating the loss as a function of space is to integrate over non-rectangular regions. For example, imagine we want to calculate the total absorption within a spherical region. In such cases, we can calculate the loss as a function of space, then integrate over the spherical volume. It is also possible to define an integration filter using the material properties to define whether a point should be included or not included in the integration. See the Advanced method - multiple materials page for more information on how to do this.
The following figure shows a cross section of an integration filter for a sphere centered at (0,0,-200 nm) and a radius of 200nm.
The final output of the script is the percentage of power absorbed within the spherical region at a wavelength of 500nm: ~1%
Calculating absorbed optical power - Higher accuracy
This page describes a more accurate spatial absorption calculation which minimizes interpolation errors that occur because of discontinuities in the electric field near interfaces. Consider using the advanced method if a significant fraction of the absorption occurs within the first mesh cell at interfaces.
Background
This technique is an extension of the simple method described on the previous page. For reasons given on the previous page, the advanced method calculates the absorption profile using the relation
Pabs=−0.5 ω|E|2 imag(ε)
The difference between the two techniques is the monitor interpolation option.
An inherent part of the FDTD method is that during the simulation, each field component (Ex, Ey, Ez) is calculated at a different location within the Yee cell. If a monitor collects raw field data from the simulation, post processing becomes complicated due to the fact that each field component is known at a different location. Even simple calculations like |E|^2 = |Ex|^2+ |Ey|^2+ |Ez|^2 are not trivial since Ex, Ey, and Ez are not known at the same location.
To simplify post-processing of simulation data, the default behavior of monitors in FDTD is to interpolate all field components back to a common set of points (the origin of the Yee cell). This makes post-processing simple, since all field components are known at a common point.
In most situations, the interpolation errors that occur when interpolating each field component back to the origin of the Yee cell is negligible. Unfortunately, the absorption calculation can be quite sensitive to these interpolation errors. The interpolation errors only become significant when:
- The electric field component normal to the structure surface is large AND
- a significant fraction of the total absorption occurs within the first set of mesh cells at the structure surface.
From Maxwell's equations, we know that the field component normal to the structure surface will be discontinuous at the interface. Also, it is obvious that interpolating near a discontinuity can lead to errors. In the absorption calculations, the interpolation error typically makes the electric field just inside the absorbing material much larger than it really is. Since absorption is proportional to electric field intensity, this will tend to overestimate the total absorption.
Note that these errors only occur for the mesh cells which cross the material interface. If the absorption in the mesh cells that cross the interface is a small fraction of the total absorption, then the overall error will be small. However, if most of the absorption is concentrated near the surface (as is the case with metals), then the overall error can be significant.
| Note: Conformal mesh The material data has an index that is on the same order of magnitude as the background index. Also, we are able to simulate this problem at a reasonably small mesh size of 5nm or less. For this reason, we can change the mesh refinement option to "conformal variant 1" to take full advantage of the conformal meshing feature, even for a metals. |
Example
The simulation file usr_absorption_advanced.fsp can be used to calculate the power absorption per unit volume of a gold sphere illuminated by a focused beam.
The total absorption is calculated in three ways:
- A box of 2D monitors around the sphere calculates the net power absorbed within the box. This is an easy, efficient and accurate technique to calculate total absorption within a rectangular region.
- The simple absorption analysis group described in the previous page. The monitors in this group use the standard field interpolation option. It will not give accurate results in this simulation because the field component normal to the surface of the gold sphere is large. Errors in the absorption data near the surface are significant since most of the absorption occurs near the surface.
- The advanced absorption analysis group. The monitors in this group do not use monitor interpolation. This gives more accurate absorption data near the sphere surface.
Run the simulation, then run the associated analysis script. It will create a figure showing the absorption spectrum as calculated by the three techniques. The analytic solution (calculated from Mie theory) is also shown.
There are a few interesting points to notice:
- The Net power into box and Advanced techniques are very similar.
- The absorption calculated by the standard technique is quite different than the result from the other techniques The standard technique is not appropriate for calculating the spatial absorption profile in metals because it is very sensitive to interpolation errors of discontinuous field components at the metallic interfaces.
- All three techniques are somewhat different than the analytic result. This example was setup to run quickly, and not to provide the most accurate result. A smaller mesh size, larger simulation region, and more PML layers would be required to get closer agreement with the theory.
| Tip: Size of the monitor To get the most accurate results, the analysis object should be large enough to completely surround the absorbing material, plus at least 1 mesh cell of non-absorbing region around the edges of the analysis object |
| Tip: Periodic and Bloch boundary conditions To get the most accurate results when using this analysis object with Periodic and Bloch boundary conditions, a small change should be made within the Analysis script. For details, see the analysis script comments or contact Lumerical technical support. |
| Note: Higher accuracy simulation The following figure shows the absorption data from a simulation that was run with higher accuracy settings. As you can see, the FDTD results are converging to the analytic solution. The simple Net power box and Advanced technique are closest to the theory.
To increase the accuracy of the example file, make the following changes:
Warning: The memory requirements and file size will be very large with these changes. |
Calculating absorbed optical power - Higher accuracy method with multiple materials
This example shows how to calculate the power absorbed in a specific material when there are two (or more) dispersive materials in a simulation and objects are of arbitrary shapes where it would be difficult to define the a spatial filter. This example applied both 2D and 3D cases.
In this example of usr_absorption_advanced_material.fsp a silver particle of arbitrary shape is embedded in a silicon substrate. Suppose we are interested in the power absorbed in the silver particle.
The advanced absorbed power monitor group returns the power absorbed as a function of space, ie Pabs(x,y,z). Once the project file is run, executing the usr_absorption_advanced_material.lsf file will multiply Pabs(x,y,z) by 1 if the (x,y,z) values are inside the silver particle and 0 if they are outside.
We can use the material properties of Silicon and Silver to tell if a specific point (x,y,z) lies inside particle or not.
The script, uses the index monitor inside the Pabs box to create the filters by comparing both the real and imaginary parts the index of every point in the simulation to those of Silver and Silicon. The figure below shows the filter generated for Silicon:
Running the script will generate the following plots for total absorption, and absorption spectra in the two materials.
This technique can be extended to 3D simulations as well.
Calculating absorbed optical power - from divergence of Poynting vector
This section describes how to calculate the absorption per unit volume from the divergence of the Poynting vector. Divergence calculations tend to be very sensitive to numerical problems, which makes this technique inferior to the other techniques described in this section.
The loss per unit volume can be calculated from the Poynting vector by
dLdV=12real(∇⋅P)
Consider the example file usr_absorption_divergence.fsp. The simulation consists of a gold nanoparticle in a background material that has an index of 1.5 + i0.05. A TFSF source is used to excite the simulation and a 3D power monitor is used to measure the Poynting vector.
We want to measure the power absorbed in the gold nanoparticle. If the background was not absorbing, a box of 2D power monitors around the particle could be used to measure the total absorbed power. In this case, since the background is also absorbing, the box technique will not work. Instead, we will measure the Poynting vector as a function of (x,y,z), and calculate the loss with the above formula. The loss per volume can be integrated over the volume of the gold sphere, giving the total power absorbed by the gold.
Load the file usr_absorption_divergence.fsp, and run the simulation. When the simulation is complete, run the script usr_absorption_divergence.lsf. It will calculate the total power absorbed within the box of monitors using two techniques. First, with a box of 2D monitors measuring the transmission through each surface. Second, by calculating the loss per unit volume, then integrating over the entire volume of the box. Both techniques give 48% absorption.
Next, the loss within the gold particle is calculated by integrating the loss over the volume of the sphere. It is quite small, only 0.1% of the total power. This seems reasonable after looking at a slice of the loss per unit volume in the XY plane. There is only a small amount of absorption in the particle very close to its surface. The matching slice of the integration filter is also shown. The power absorbed in any arbitrary volume can be calculated by modifying the integration filter.
| Note: 3D monitors can require large amounts of memory and computation time. Therefore, it is best to keep them as small as possible. Similarly, keep the number of frequency points small, and only record the field components that are necessary. In this example, only the Poynting vector was saved, not the E and H field components. |
Determining the direction of power flow
In some cases, it is interesting to understand if power is flowing forwards or backwards in a particular plane of a simulation. In this example, we show how this information can be obtained for a simulation with periodic or Bloch boundary conditions when there is a planar region of homogeneous, low-loss dielectric material.
This information cannot easily be determined from considering the Poynting vector alone. If we consider a simple example of a plane wave incident upon a perfectly reflecting boundary, the Poynting vector and the total transmitted power through a plane will be zero because the total normalized transmitted power is 1 (incident wave) - 1 (reflected wave) = 0. This is not easy to distinguish this from the case of zero fields everywhere, which also gives a total normalized transmitted power of 0.
In this example, we use the power direction analysis object which uses farfield and grating projections to decompose the measured E and H fields onto a basis state of forward and backward propagating plane waves. The forward power is then calculated by summing the power in all the forward propagating plane waves, and similarly for the backward propagating power. This method has some numerical errors that are largest when we have plane waves propagating at very steep angles to the measurement plane. For this reason, the forward transmission minus the backwards transmission may be different from the net transmission, measured simply by integrating the Poynting vector, by a few percent. The power direction object can be inserted from the object library in the optical power section.
Example
In this example, we consider light incident in a medium of dielectric index 1.4. The light propagates towards a gold ellipsoid on a substrate of index 2. Periodic boundary conditions are used, meaning that we are considering the light incident on a periodic array of these gold nanoparticles. The period is 1 micron in both x and y. We place a special analysis group for this calculation between the source and the particle so that it measures both the incident and reflected field. For the purpose of testing, we also place a reflection monitor behind the source. The analysis group contains a planar power monitor and also a small index monitor because it is necessary to know the refractive index where the fields were measured. It is important to note that this calculation is only valid if:
- The imaginary part of the refractive index is very small or zero. We ignore the imaginary part in this calculation.
- The material is completely uniform over the plane where the monitor is recording the field.
- Periodic or Bloch boundary conditions are used. This method could be extended to cases where PML is used as long as the fields go to zero near the edges of the simulation. Please contact Lumerical support if you would like to explore this case.
After running the simulation usr_power_direction.fsp , please run the script usr_power_direction_test.lsf . This analysis can take about a minute because there are a large number of far field projections involved and we record data at several different wavelengths - a new farfield projection is required for each wavelength. The script will create two figures. The first figure shows the reflected power measured using the power direction group, which is placed between the source and the nanoparticle. We see that we have successfully calculated the amount of backward propagating power ("new method"), which is close the reflected power measured in the usual way ("standard method"), although there is a small numerical error. The second figure compares the total transmission through the monitor measured in 3 ways:
- by using T_total, calculated from the transmission function
- by calculating R-1 where R is the reflection calculated in the usual way
- by taking the difference of T_positive and T_negative
We see that all methods agree reasonably well, but there is the largest error in T_positive-T_negative because this method involves the plane wave decomposition that introduces some more numerical error.
Calculating the net power flow with a Power transmission box
This page provides a simple analysis group that calculates the net power flow out of a rectangular volume within a simulation. The files in this section were created using FDTD, but the same analysis group can be found in the Component library in MODE.
Lumerical provides many built in analysis groups in our object library. Please press this button
to open the online library of analysis groups and select optical power category to insert transmission boxes.
The monitor group in the associated files (usr_transmission_3Dbox.fsp and usr_transmission_2Dbox.fsp ) calculates the net power flow out of the box of monitors. To use this group to calculate net power flow into the box, just multiply the result by -1.
For example, in the associated files, we have a small silicon particle in a focused beam. The box of monitors can be used to measure the power absorbed by the particle. After running the simulation, use the following commands to plot the absorption vs wavelength. Once the visualizer is open, you may want to use the "-Real" option to make the numbers positive, since power is flowing into the box.
runanalysis;
net_power = getresult("trans_box","T");
visualize(net_power);
| Note: Symmetry These monitor boxes work properly when symmetry boundary conditions are used. When symmetry boundary conditions are used, the monitor group assumes there will be equal power flows on both sides of the plane of symmetry. |
| Note: Field data not saved The monitors in these groups are set to record the net power only, to minimize the amount of data saved to file. |
| Note: Spherical transmission monitors The files usr_transmission_3Dsphere.fsp and usr_transmission_2Dsphere.fsp contain spherical transmission monitor objects. These are very similar to the transmission boxes described above, but they use a series of point monitors to create a spherical surface, rather than using 6 2D monitors to create a box. Due to the additional interpolation required by these monitors, they tend to be less accurate than box monitors. They are also not designed to work with symmetry boundary conditions and can't be used to do far field projections. For these reasons, we recommend using box monitors rather than spherical monitors.
|
Calculating the effective mode area of a waveguide mode
This page provides a simple analysis group that calculates the effective mode area.
The effective mode area, A, is the ratio of a mode's total energy density per unit length and its peak energy density
A=1max{W(r)}∫A∞W(r)dA
where W(r) is the energy density,
W(r)=12Re{d[ωε(r)]dω}|E(r)|2+12μ0|H(r)|2
Example
The simulation file usr_effective_mode_area.fsp contains a simple silicon on insulator (SOI) waveguide shown in the image above. To get the effective mode area for the injected mode, first run the simulation. Then edit the effective mode area analysis group and press the RUN ANALYSIS button in the ANALYSIS-->SCRIPT tab. As can be seen in the image below, the analysis script output will contain the calculated effective mode area.
In the screen shot above, you can see that the analysis script contained in the usr_effective_mode_area.fsp simulation uses two methods to compute the effective mode area. The two methods only differ in how the derivative in the energy density is computed. Method 1 uses central differences on the derivative on the left hand side of the equation below. And method 2 uses central differences on the derivative on the right hand side of the equation.
Re{d[ωε(r)]dω}=Re{ωd[ε(r)]dω+ε(r)}
Calculating the modal volume of a cavity mode
This page provides simple example using the modal volume analysis object that can be used to calculate modal volume of a confined mode. This object is available in the object library in the resonators section. A modal area analysis object is also available from the object library for 2D simulations.
Modal volume calculation
This analysis object provides three simple ways to calculate the modal volume, defined below.
| Calculation Type | Equation |
|---|---|
| 1 | V=∫(E2>12E2)dV |
| 2 | V=∫εE2dVmax(εE2) |
| 3 (default) | V=(∫H2dV)2∫H4dV |
It is worth noting that more sophisticated definitions of modal volume exist. If these definitions are not sufficient for your needs, these examples can be used as a starting point for you to create your own custom modal volume analysis object.
Other analysis
Impulse Analysis - Calculating the response to an arbitrary time signal using the impulse response
This topic explains how to calculate the time domain response of a system to an arbitrary source time signal. The response is calculated by using the impulse response of the system (as calculated from a simulation that uses the standard source pulse).
This type of post-processing analysis is often more efficient than directly running an FDTD simulation with the arbitrary source pulse. For example, once the impulse response is known, it is possible to calculate the response to many different inputs without needing to run any additional simulations.
Methodology
The impulse response of the system is obtained by running a 'standard' simulation using the default source pulse. Important points to consider include:
- The source must be setup to excite the system over the full bandwidth of interest (ie. it must cover the full spectrum of the desired custom source signal).
- A frequency-domain monitor must be setup to measure the impulse response of the system over the same range of frequencies. It is important to understand that the frequency monitors record the impulse response of the system. No additional processing of the data is required to obtain the impulse response.
- As with any simulation, it is important that the simulation time be sufficiently long so that the time domain fields decay to ~0 by the end of the simulation.
- Note that the monitors in this example file have been setup to always record data from 0.7-1.3um, regardless of the source wavelength range.
Once the impulse response has been calculated by a frequency monitor, it is possible to calculate the system's response to an arbitrary input based on some simple fourier analysis.
- The user must specify the desired input signal as a function of time: Signal(t)
- Calculate the input signal's frequency domain spectrum with a fourier transform: Signal(f)
- Multiply the signal spectrum and the impulse response: Impulse(f) * Signal(f)
- Inverse transform the spectrum*impulse data to obtain the systems response to the desired input: Response(t).
Example
To reproduce the following examples, open the associated simulation file (.fsp) and run the script file (.lsf). The goal of this example is to determine the response of the system to a 50fs Gaussian pulse propagating through a n=3 dielectric slab using the technique described above. To verify the result, the script runs another FDTD simulation using a 50 Gaussian pulse so the calculated response can be compared to the response from the direct FDTD simulation.
| Impulse response of system for fields transmitted through the structure from 0.7-1.3um. Oscillation is due to constructive/destructive interference as fields pass through the dielectric slab. |
|
| Spectrum of 50fs gaussian pulse. Notice that spectrum of pulse is fully contained within the range of known impulse response data (0.7-1.3um). |
|
| Time domain response of system Input signal: Desired input signal. This is typically specified by the user. Response - from simulation: The response of the system as calculated by running a full FDTD simulation. For verification purposes. Response - from impulse analysis: The response of the system as calculated from the source spectrum and impulse response data. Notice that this data is very similar to the response as calculated from the verification simulation (green line), but did not require an additional FDTD simulation. |
|
6260
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
