Finite Difference Time Domain Solver
=====================================
.. contents::
   :local:
   :depth: 2
   :backlinks: top
   
The Finite-Difference Time-Domain (FDTD) solver is a powerful electromagnetic simulation tool for analyzing integrated photonic devices. This guide will walk you through setting up and running your first FDTD simulation using a straight waveguide example.

Basic Workflow
--------------

Every FDTD simulation follows these steps:

1. **Define materials** with their optical properties
2. **Create layer stack** describing the vertical structure
3. **Set up device geometry** from a GDS layout file or function
4. **Configure ports** for mode injection and monitoring
5. **Add monitors** to capture field data or S-parameters
6. **Set simulation parameters** including time, mesh, and boundaries
7. **Run simulation** and analyze results

Example: Straight Waveguide Simulation
---------------------------------------

This example demonstrates the essential steps to simulate a straight silicon waveguide and extract S-parameters.

Step 1: Import Required Modules
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

.. code-block:: python

   from pyOptiShared.DeviceGeometry import DeviceGeometry
   from pyOptiShared.LayerInfo import LayerStack
   from pyOptiShared.Material import ConstMaterial
   from pyFDTDKernel.pyFDTDSolver import pyFDTDSolver
   import gdstk

Step 2: Define Materials
^^^^^^^^^^^^^^^^^^^^^^^^^

Create materials with constant refractive indices:

.. code-block:: python

   # Silicon dioxide (n=1.444) for substrate and cladding
   sio2_mat = ConstMaterial(mat_name="SiO2", epsReal=1.444**2, color='lightgreen')
   
   # Silicon (n=3.48) for waveguide core
   si_mat = ConstMaterial(mat_name="Si", epsReal=3.48**2, color='lightblue')
   
   # Air (n=1.0) for background
   air_mat = ConstMaterial(mat_name="Air", epsReal=1**2, color='lightyellow')

**Note**: The dielectric constant (epsilon) is the square of the refractive index: ε = n\ :sup:`2`

Step 3: Build the Layer Stack
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

The layer stack defines your device's vertical structure:

.. code-block:: python

   layer_stack = LayerStack()
   
   # Add silicon layer: 220nm thick starting at z=0
   layer_stack.addLayer(
       name="L1", 
       number=1, 
       thickness=0.22, 
       zmin=0.0,
       material=si_mat, 
       cladding=air_mat,
       sideWallAng=0
   )
   
   # Set background (above structure) and substrate (below structure)
   layer_stack.setBGandSub(background=air_mat, substrate=sio2_mat)

Step 4: Create Device Geometry
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

You can load geometry from a GDS file or define it programmatically. Here we'll use gdstk to generate "waveguide.gds" file which will be used to define the geometry:

.. code-block:: python

   length = 10
   width = 0.5
   layer_core = 1

   output_filename = "waveguide.gds"
   lib = gdstk.Library()

   strt_wg = lib.new_cell("Straight")
   vertices = [(0, -width/2), (length, -width/2), (length, width/2), (0, width/2)]
   strt_wg.add(gdstk.Polygon(vertices, layer=layer_core))

   lib.write_gds(output_filename)

   device_geometry = DeviceGeometry()
   device_geometry.SetFromGDS(
       layer_stack=layer_stack,
       gds_file='waveguide.gds',
       buffers={'x': 1.5, 'y': 1.5, 'z': 1.5}  # 1.5 μm padding on all sides
   )
   
   # Automatically detect ports along the x-direction
   device_geometry.SetAutoPortSettings(
       direction='x',
       port_buffer=1.0
   )
   
   device_geometry.PlotGDS()
   device_geometry.PlotGDSXSection(direction='y',pos=2)


**Port Detection**: The solver automatically identifies waveguide ports at domain boundaries. It is also possible to define ports manually if auto port detection fails.

.. grid:: 1

    .. grid-item-card::  GDS Mask

        .. figure:: ../_static/FDTD_GettingStarted_GDS.png
           :width: 100%
           :align: center

    .. grid-item-card::  Waveguide Cross Section

        .. figure:: ../_static/FDTD_GettingStarted_GDS_xsection.png
           :width: 100%
           :align: center


Step 5: Configure the FDTD Solver
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

Set up the solver with excitation and monitoring:

.. code-block:: python

   # Define wavelength range
   lmin = 1.5   # Minimum wavelength (μm)
   lmax = 1.6   # Maximum wavelength (μm)
   lcen = (lmax + lmin) / 2  # Center wavelength
   npts = 21    # Number of frequency points
   
   fdtd_solver = pyFDTDSolver()
   
   # Configure port excitation with Gaussian pulse
   fdtd_solver.SetPorts(
       profile="gaussian-pw",      # Gaussian pulse
       lcenter=lcen,
       lmin=lmin,
       lmax=lmax,
       npts=npts,
       mode_indices=0,             # Excite fundamental mode (TE0)
       symmetries='1x1'            # Use symmetry for single port
   )

**Key Parameters**:

- **profile**: Type of excitation (gaussian-pw for broadband analysis)
- **mode_indices**: Which mode to excite (0 = fundamental mode)
- **symmetries**: Exploit device symmetry to reduce simulation time

Step 6: Add Monitors
^^^^^^^^^^^^^^^^^^^^

Add a DFT (Discrete Fourier Transform) monitor to capture frequency-domain fields:

.. code-block:: python

   fdtd_solver.AddDFTMonitor(
       mon_type="2d-z-normal",     # 2D monitor in XY plane
       z0=0.11,                     # Z-position (center of waveguide)
       name="MyDFTMonitor1",
       lmin=lmin,
       lmax=lmax,
       npts=npts,
       save_ex=True,                # Save Ex field component
       save_ey=True,                # Save Ey field component
       save_hz=True                 # Save Hz field component
   )

**Monitor Types**:

- **2d-z-normal**: XY plane (for waveguides propagating in X or Y)
- **2d-x-normal**: YZ plane
- **2d-y-normal**: XZ plane

Step 7: Set Simulation Parameters
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

Configure mesh resolution, simulation time, and output:

.. code-block:: python

   fdtd_solver.SetSimSettings(
       sim_time=350,                    # Total simulation time (ps)
       space_step=0.05,                 # Mesh size: 50nm
       subpixel_level=2,                # Subpixel averaging level
       save_path="./results",
       results_filename='waveguide',
       device_geometry=device_geometry,
       auto_shutoff_limit=1e-2,         # Auto-stop when fields decay
       export_mat_grid=True             # Export permittivity grid
   )

**Critical Parameters**:

- **space_step**: Mesh resolution (smaller = more accurate but slower)
- **subpixel_level**: Material averaging (1=off, 2=recommended)
- **sim_time**: Duration in time units (depends on device size)
- **auto_shutoff_limit**: Stops simulation when field energy falls below threshold

Step 8: Run Simulation
^^^^^^^^^^^^^^^^^^^^^^

.. code-block:: python

   # Run the simulation
   results = fdtd_solver.Run()
   
   # Plot S-parameters
   results.PlotSParameters(snp='ALL', plot_type='power')
   results.PlotSParameters(snp='ALL', plot_type='phase')


.. grid:: 1

    .. grid-item-card::  S-Parameters (Power)

        .. figure:: ../_static/FDTD_GettingStarted_sparam_pow.png
           :width: 100%
           :align: center

    .. grid-item-card::  S-Parameters (Phase)

        .. figure:: ../_static/FDTD_GettingStarted_sparam_phi.png
           :width: 100%
           :align: center

Complete Example Code
^^^^^^^^^^^^^^^^^^^^^

Here is the complete code you can copy and run (:download:`FDTD_getting_started.py <../_static/examples/GettingStarted/FDTD_getting_started.py>`):

.. literalinclude:: ../_static/examples/GettingStarted/FDTD_getting_started.py
   :language: python

Understanding the Results
-------------------------

The FDTD solver returns an ``FDTDSimResults`` object with several useful components:

S-Parameters
^^^^^^^^^^^^

S-parameters quantify how much power is transmitted and reflected:

- **S11**: Reflection at Port 1 (how much power bounces back)
- **S21**: Transmission from Port 1 to Port 2 (how much power passes through)
- Power transmission: \|S21\|\ :sup:`2`
- Power reflection: \|S11\|\ :sup:`2`

Field Data
^^^^^^^^^^

DFT monitors store frequency-domain field components:

- Access fields: ``results.runs[0].dftmonitors["MyDFTMonitor1"].Get('Hz')``
- Plot field profiles to visualize mode propagation

Results Storage
^^^^^^^^^^^^^^^

All data is saved in HDF5 format (``results/waveguide.hdf5``):

.. code-block:: python

   from pyOptiShared.SimResults import FDTDSimResults
   
   # Load results from file
   results = FDTDSimResults()
   results.loadHDF5('results/waveguide.hdf5')

Tips for Success
----------------

1. **Mesh Resolution**: Start with 50nm (0.05 μm) and refine if needed. Smaller mesh = more accuracy but longer runtime.

2. **Subpixel Averaging**: Always use ``subpixel_level=2`` for curved or diagonal features. This dramatically improves accuracy without increasing mesh density.

3. **Simulation Time**: Ensure fields have time to propagate through the device and decay. Monitor the auto-shutoff to verify convergence.

4. **Boundary Padding**: Use at least 1 μm buffers around devices to prevent reflections from PML boundaries.

5. **Mode Selection**: ``mode_indices=0`` excites the fundamental mode. For multimode devices, use higher indices (1, 2, etc.).

6. **Wavelength Range**: Choose range based on your application. More frequency points give smoother spectra but take longer.

Advanced Features
-----------------

Port Symmetries
^^^^^^^^^^^^^^^

For symmetric devices (crossings, splitters), exploit symmetry to simulate only one port:

.. code-block:: python

   fdtd_solver.SetPorts(
       ...,
       symmetries='4x'  # Predefined 4-port crossing symmetry
   )

Time Monitors
^^^^^^^^^^^^^

Capture time-domain field evolution for animations:

.. code-block:: python

   fdtd_solver.AddTimeMonitor(
       mon_type="2d-z-normal",
       z0=0.11,
       name="MyTimeMonitor1",
       npts=200  # Number of time snapshots
   )

High Order Modes
^^^^^^^^^^^^^^^^^^

Excite specific modes in different ports:

.. code-block:: python

   fdtd_solver.SetPorts(
       ...,
       mode_indices=[0, 2]  # Port 1: mode 0, Port 2: mode 2
   )

Next Steps
----------

Now that you understand the basics, you can:

- Simulate bends, splitters, and other photonic components
- Analyze wavelength-dependent transmission and reflection
- Export field profiles for visualization
- Optimize device geometries for target performance

For advanced examples including 90-degree bends, polarization converters, and waveguide crossings, consult the full documentation :ref:`fdtd_api`.