# Germanium resonator benchmark#

In this example, we reproduce the findings of Campione et al. (2016), which is linked here.

This notebook was originally developed and written by Romil Audhkhasi (USC).

The paper investigates the resonances of Germanium structures by measuring their transmission spectrum under varying geometric parameters.

The paper uses a different commercial finite-difference time-domain software, which matches the result from Tidy3D.

To do this calculation, we use a broadband pulse and frequency monitor to measure the flux on the opposite side of the structure.

[1]:

# standard python imports
import numpy as np
import matplotlib.pyplot as plt

# tidy3D import
import tidy3d as td
from tidy3d import web


[10:32:33] WARNING  This version of Tidy3D was pip installed from the         __init__.py:100
'tidy3d-beta' repository on PyPI. Future releases will be

           INFO     Using client version: 1.9.0                               __init__.py:115


## Set Up Simulation#

[2]:

Nfreq = 1000
wavelengths = np.linspace(8, 12, Nfreq)
freqs = td.constants.C_0 / wavelengths
freq0 = freqs[len(freqs) // 2]
freqw = freqs[0] - freqs[-1]

# Define material properties
n_BaF2 = 1.45
n_Ge = 4
BaF2 = td.Medium(permittivity=n_BaF2**2)
Ge = td.Medium(permittivity=n_Ge**2)


[3]:

# space between resonators and source
spc = 8

# geometric parameters
Px = Py = P = 4.2
h = 2.53
L1 = 3.036
L2 = 2.024
w1 = w2 = w = 1.265

# resolution (should be commensurate with periodicity)
dl = P / 32


[4]:

# total size in z and [x,y,z]
Lz = spc + h + h + spc
sim_size = [Px, Py, Lz]

# BaF2 substrate
substrate = td.Structure(
geometry=td.Box(
center=[0, 0, -Lz / 2],
size=[td.inf, td.inf, 2 * (spc + h)],
),
medium=BaF2,
name="substrate",
)

# Define structure

cell1 = td.Structure(
geometry=td.Box(
center=[(L1 / 2) - L2, -w1 / 2, h / 2],
size=[L1, w1, h],
),
medium=Ge,
name="cell1",
)

cell2 = td.Structure(
geometry=td.Box(
center=[-L2 / 2, w2 / 2, h / 2],
size=[L2, w2, h],
),
medium=Ge,
name="cell2",
)


[5]:

# time dependence of source
gaussian = td.GaussianPulse(freq0=freq0, fwidth=freqw)

# plane wave source
source = td.PlaneWave(
source_time=gaussian,
size=(td.inf, td.inf, 0),
center=(0, 0, Lz / 2 - spc + 2 * dl),
direction="-",
pol_angle=0,
)

# Simulation run time.  Note you need to run a long time to calculate high Q resonances.
run_time = 3e-11


[6]:

# monitor fields on other side of structure (substrate side) at range of frequencies
monitor = td.FluxMonitor(
center=[0.0, 0.0, -Lz / 2 + spc - 2 * dl],
size=[td.inf, td.inf, 0],
freqs=freqs,
name="flux",
)



## Define Case Studies#

Here we define the two simulations to run

• With no resonator (normalization)

• With Ge resonator

[7]:

grid_spec = td.GridSpec(
grid_x=td.UniformGrid(dl=dl),
grid_y=td.UniformGrid(dl=dl),
grid_z=td.AutoGrid(min_steps_per_wvl=32),
)

# normalizing run (no Ge) to get baseline transmission vs freq
# can be run for shorter time as there are no resonances
sim_empty = td.Simulation(
size=sim_size,
grid_spec=grid_spec,
structures=[substrate],
sources=[source],
monitors=[monitor],
run_time=run_time / 10,
boundary_spec=td.BoundarySpec(x=td.Boundary.periodic(), y=td.Boundary.periodic(), z=td.Boundary.pml()),
)

# run with Ge nanorod
sim_actual = td.Simulation(
size=sim_size,
grid_spec=td.GridSpec.uniform(dl=dl),
structures=[substrate, cell1, cell2],
sources=[source],
monitors=[monitor],
run_time=run_time,
boundary_spec=td.BoundarySpec(x=td.Boundary.periodic(), y=td.Boundary.periodic(), z=td.Boundary.pml()),
)


[8]:

# Structure visualization in various planes

fig, (ax1, ax2, ax3) = plt.subplots(1, 3, figsize=(10, 6))
sim_actual.plot_eps(x=0, ax=ax1)
sim_actual.plot_eps(y=-0.1, ax=ax2)
sim_actual.plot_eps(z=0.1, ax=ax3)
plt.show()



## Run Simulations#

[9]:

# run all simulations, take about 2-3 minutes each with some download time
batch = web.Batch(simulations={"norm": sim_empty, "actual": sim_actual})
batch_data = batch.run(path_dir="data")


[10:32:34] INFO     Auto meshing using wavelength 10.0020 defined from       grid_spec.py:510
sources.

           INFO     Using Tidy3D credentials from stored file.                     auth.py:77

[10:32:35] INFO     Authentication successful.                                     auth.py:37

           INFO     Created task 'norm' with task_id                            webapi.py:131
'62f6a87e-cb5f-48c4-a253-30fd578b10fb'.

[10:32:36] INFO     Created task 'actual' with task_id                          webapi.py:131
'cedbcea8-d28c-4146-a2c1-fdd13019e286'.

[10:32:37] Started working on Batch.                                         container.py:383

[10:33:13] Batch complete.                                                   container.py:417


The normalizing run computes the transmitted flux for an air -> SiO2 interface, which is just below unity due to some reflection.

While not technically necessary for this example, since this transmission can be computed analytically, it is often a good idea to run a normalizing run so you can accurately measure the change in output when the structure is added. For example, for multilayer structures, the normalizing run displays frequency dependence, which would make it prudent to include in the calculation.

[10]:

batch_data = batch.load(path_dir="data")
transmission = batch_data["actual"]["flux"].flux / batch_data["norm"]["flux"].flux
reflection = 1 - transmission


[10:33:14] INFO     downloading file "output/monitor_data.hdf5" to              webapi.py:673
"data/cedbcea8-d28c-4146-a2c1-fdd13019e286.hdf5"

[10:33:15] INFO     loading SimulationData from                                 webapi.py:472
data/cedbcea8-d28c-4146-a2c1-fdd13019e286.hdf5

           WARNING  Simulation final field decay value of 4.07e-05 is greater   webapi.py:479
than the simulation shutoff threshold of 1e-05. Consider
simulation again with large run_time duration for more
accurate results.

           INFO     downloading file "output/monitor_data.hdf5" to              webapi.py:673
"data/62f6a87e-cb5f-48c4-a253-30fd578b10fb.hdf5"

           INFO     loading SimulationData from                                 webapi.py:472
data/62f6a87e-cb5f-48c4-a253-30fd578b10fb.hdf5

[11]:

# plot transmission, compare to paper results, look similar
fig, ax = plt.subplots(1, 1, figsize=(6, 4.5))
plt.plot(wavelengths, reflection, "k", label="R")
plt.plot(wavelengths, transmission, "r--", label="T")
plt.xlabel("wavelength ($\mu m$)")
plt.ylabel("Magnitude")
plt.xlim([8.8, 12])
plt.ylim([0.0, 1.0])
plt.legend()
plt.show()


[ ]: