Further Analogue-Enhanced Cosimulation including SPICE
#
This example demonstrates the modelling of multi-physical component interconnection and system design.
[ ]:
from gdsfactory.components import mzi2x2_2x2_phase_shifter
import hdl21 as h
import pandas as pd
import numpy as np
import piel
import sax
import sys
Start from gdsfactory
#
Let us begin from the mzi2x2_2x2_phase_shifter
circuit we have been analysing so far in our examples, but this time, we will explore further aspects of the electrical interconnect, and try to extract or extend the simulation functionality. Note that this flattened netlist is the lowest level connectivity analysis, for very computationally complex extraction. Explore the electrical circuit connectivity using the keys below.
Note that this netlist just gives us electrical ports and connectivity for this component.
mzi2x2_2x2_phase_shifter_netlist_flat["connections"]The most basic model of this phase-shifter device is treating the heating element as a resistive element. We want to explore how this resistor model affects the rest of the photonic design circuitry. You can already note something important if you explored the netlist with sufficient care and followed the previous examples regarding the straight_x_top
or sxt
heater component.
[ ]:
mzi2x2_2x2_phase_shifter_netlist = mzi2x2_2x2_phase_shifter().get_netlist(
exclude_port_types="optical"
)
mzi2x2_2x2_phase_shifter_netlist["instances"]["sxt"]
{'component': 'straight_heater_metal_undercut',
'info': {'resistance': None},
'settings': {'cross_section': 'strip',
'length': 200.0,
'with_undercut': False}}
So this top heater instance info
instance definition, it already includes a resistance
field. However, in the default component definition, it is defined as None
. Let us give some more details about our circuit, and this would normally be provided by the PDK information of your foundry.
[ ]:
from piel import straight_heater_metal_simple
import functools
# Defines the resistance parameters
our_resistive_heater = functools.partial(
straight_heater_metal_simple, ohms_per_square=2
)
our_resistive_mzi_2x2_2x2_phase_shifter = mzi2x2_2x2_phase_shifter(
straight_x_top=our_resistive_heater,
)
our_resistive_mzi_2x2_2x2_phase_shifter.plot_widget()
[ ]:
our_resistive_mzi_2x2_2x2_phase_shifter_netlist = (
our_resistive_mzi_2x2_2x2_phase_shifter.get_netlist(exclude_port_types="optical")
)
our_resistive_mzi_2x2_2x2_phase_shifter_netlist["instances"]["sxt"]
{'component': 'straight_heater_metal_undercut',
'info': {'resistance': 1000.0},
'settings': {'cross_section': 'strip',
'length': 200.0,
'ohms_per_square': 2,
'with_undercut': False}}
Now we have a resistance parameter directly connected to the geometry of our phase shifter topology, which is an incredibly powerful tool. In your models, you can compute or extract resisivity and performance data parameters of your components in order to create SPICE models from them. This is related to the way that you perform your component design. piel
can then enable you to understand how this component affects the total system performance.
We can make another variation of our phase shifter to explore physical differences.
[ ]:
our_short_resistive_mzi_2x2_2x2_phase_shifter = mzi2x2_2x2_phase_shifter(
straight_x_top=our_resistive_heater,
length_x=100,
)
our_short_resistive_mzi_2x2_2x2_phase_shifter.show()
our_short_resistive_mzi_2x2_2x2_phase_shifter.plot_widget()
You can also find out more information of this component through:
[ ]:
our_short_resistive_mzi_2x2_2x2_phase_shifter.named_references["sxt"].info
{‘resistance’: 500.0}
So this is very cool, we have our device model giving us electrical data when connected to the geometrical design parameters. What effect does half that resistance have on the driver though? We need to first create a SPICE model of our circuit. One of the main complexities now is that we need to create a mapping between our component models and hdl21
which is dependent on our device model extraction. Another functionality we might desire is to validate physical electrical connectivity by
simulating the circuit accordingly.
Extracting the SPICE circuit and assigning model parameters#
We will exemplify how piel
microservices enable the extraction and configuration of the SPICE circuit. This is done by implementing a SPICE netlist construction backend to the circuit composition functions in sax
, and is composed in a way that is then integrated into hdl21
or any SPICE-based solver through the VLSIR
Netlist
.
The way this is achieved is by extracting all the instances
, connections
, ports
, models
which is essential to compose our circuit using our piel
SPICE backend solver. It follows a very similar principle to all the other sax
based circuit composition tools, which are very good.
[ ]:
our_resistive_heater_netlist = our_resistive_heater().get_netlist(
allow_multiple=True, exclude_port_types="optical"
)
# our_resistive_mzi_2x2_2x2_phase_shifter_netlist = our_resistive_mzi_2x2_2x2_phase_shifter.get_netlist(exclude_port_types="optical")
# our_resistive_heater_netlist
We might want to extract our connections of our gdsfactory netlist, and convert it to names that can directly form part of a SPICE netlist. However, to do this we need to assign what type of element each component each gdsfactory instance is. We show an algorithm that does the following in order to construct our SPICE netlist:
Parse the gdsfactory netlist, assign the electrical ports for the model. Extract all instances and required models from the netlist.
Verify that the models have been provided. Each model describes the type of component this is, how many ports it requires and so on. Create a
hdl21
top level module for every gdsfactory netlist, this is reasonable as it is composed, and not a generator class. This generates a large amount of instantiatedhdl21
modules that are generated fromgenerators
.Map the connections to each instance port as part of the instance dictionary. This parses the connectivity in the
gdsfactory
netlist and connects the ports accordingly.
piel
does this for you already:
[ ]:
our_resistive_heater_spice_netlist = piel.gdsfactory_netlist_with_hdl21_generators(
our_resistive_heater_netlist
)
our_resistive_heater_spice_netlist
This will allow us to create our SPICE connectivity accordingly because it is in a suitable netlist format using hdl21
. Each of these components in this netlist is some form of an electrical model or component. We start off from our instance definitions. They are in this format:
[ ]:
our_resistive_heater_netlist["instances"]["straight_1"]
{'component': 'straight',
'info': {'length': 320.0,
'width': 0.5,
'cross_section': 'strip_heater_metal',
'settings': {'width': 0.5,
'layer': 'WG',
'heater_width': 2.5,
'layer_heater': 'HEATER'},
'function_name': 'strip_heater_metal'},
'settings': {'cross_section': 'strip_heater_metal',
'heater_width': 2.5,
'length': 320.0},
'hdl21_model': Generator(name=Straight)}
We can compose our SPICE using hdl21
using the models we have provided. The final circuit can be extracted accordingly:
[ ]:
our_resistive_heater_circuit = piel.construct_hdl21_module(
spice_netlist=our_resistive_heater_spice_netlist
)
our_resistive_heater_circuit.instances
{'straight_1': Instance(name=straight_1 of=GeneratorCall(gen=straight)),
'taper_1': Instance(name=taper_1 of=GeneratorCall(gen=taper)),
'taper_2': Instance(name=taper_2 of=GeneratorCall(gen=taper)),
'via_stack_1': Instance(name=via_stack_1 of=GeneratorCall(gen=via_stack)),
'via_stack_2': Instance(name=via_stack_2 of=GeneratorCall(gen=via_stack))}
Note that each component is mapped into hdl21
according to the same structure and names as in the gdsfactory
netlist, if you have defined your generator components correctly. Note that the unconnected ports need to be exposed for proper SPICE composition.
[ ]:
our_resistive_heater_circuit.ports
{'e1': Signal(name='e1', width=1, desc=None),
'e2': Signal(name='e2', width=1, desc=None),
'via_stack_1__e1': Signal(name='via_stack_1__e1', width=1, desc=None),
'via_stack_1__e2': Signal(name='via_stack_1__e2', width=1, desc=None),
'via_stack_1__e4': Signal(name='via_stack_1__e4', width=1, desc=None),
'via_stack_2__e2': Signal(name='via_stack_2__e2', width=1, desc=None),
'via_stack_2__e3': Signal(name='via_stack_2__e3', width=1, desc=None),
'via_stack_2__e4': Signal(name='via_stack_2__e4', width=1, desc=None)}
Same for the signals
[ ]:
our_resistive_heater_circuit.signals
{'via_stack_1_e3': Signal(name='via_stack_1_e3', width=1, desc=None),
'via_stack_2_e1': Signal(name='via_stack_2_e1', width=1, desc=None)}
We can explore the models we have provided too:
We can also extract information for our subinstances. We can also get the netlist of each subinstance using the inherent hdl21
functionality:
[ ]:
h.netlist(
our_resistive_heater_circuit.instances["straight_1"].of, sys.stdout, fmt="spice"
)
* Anonymous `circuit.Package`
* Generated by `vlsirtools.SpiceNetlister`
*
.SUBCKT Straight__
+ e1 e2
* No parameters
rr1
+ e1 e2
+ 1000
* No parameters
.ENDS
One of the main complexities of this translation is that for the SPICE to be generated, the network has to be valid. Sometimes, in direct gdsfactory
components, there could be incomplete or undeclared port networks. This means that for the SPICE to be generated, we have to fix the connectivity in some form, and means that there might not be a direct translation from gdsfactory
. This is inherently related to the way that the construction of the netlists are generated. This means that to
some form, we need to connect the unconnected ports in order for the full netlist to be generated. The complexity is in components such as the via where there are four ports to it on each side. SPICE would treat them as four different inputs that need to be connected.
Now let’s extract the SPICE for our heater circuit:
[ ]:
h.netlist(our_resistive_heater_circuit, sys.stdout, fmt="spice")
* Anonymous `circuit.Package`
* Generated by `vlsirtools.SpiceNetlister`
*
.SUBCKT Straight__
+ e1 e2
* No parameters
rr1
+ e1 e2
+ 1000
* No parameters
.ENDS
.SUBCKT Taper__
+ e1 e2
* No parameters
rr1
+ e1 e2
+ 1000
* No parameters
.ENDS
.SUBCKT ViaStack__
+ e1 e2 e3 e4
* No parameters
rr1
+ e1 e2
+ 1000
* No parameters
rr2
+ e2 e3
+ 1000
* No parameters
rr3
+ e3 e4
+ 1000
* No parameters
rr4
+ e4 e1
+ 1000
* No parameters
.ENDS
.SUBCKT straight_heater_metal_s_b8a2a400
+ e1 e2 via_stack_1__e1 via_stack_1__e2 via_stack_1__e4 via_stack_2__e2 via_stack_2__e3 via_stack_2__e4
* No parameters
xstraight_1
+ e1 e2
+ Straight__
* No parameters
xtaper_1
+ via_stack_1_e3 e1
+ Taper__
* No parameters
xtaper_2
+ via_stack_2_e1 e2
+ Taper__
* No parameters
xvia_stack_1
+ via_stack_1__e1 via_stack_1__e2 via_stack_1_e3 via_stack_1__e4
+ ViaStack__
* No parameters
xvia_stack_2
+ via_stack_2_e1 via_stack_2__e2 via_stack_2__e3 via_stack_2__e4
+ ViaStack__
* No parameters
.ENDS
So this seems equivalent to the gdsfactory component representation. We can now continue to implement our SPICE simulation.
SPICE
Integration#
We have seen in the previous example how to integrate digital-driven data with photonic circuit steady-state simulations. However, this is making a big assumption: whenever digital codes are applied to photonic components, the photonic component responds immediately. We need to account for the electrical load physics in order to perform more accurate simulation models of our systems.
piel
provides a set of basic models for common photonic loads that integrates closely with gdsfactory
. This will enable the multi-domain co-design we like when integrated with all the open-source tools we have previously exemplified. You can find the list and definition of the provided models in:
Currently, these models do not physically represent electrically the components yet. We will do this later.
[ ]:
piel_hdl21_models = piel.models.physical.electronic.get_default_models()
piel_hdl21_models
We can extract the SPICE out of each model
[ ]:
example_straight_resistor = piel_hdl21_models["straight"]()
h.netlist(example_straight_resistor, sys.stdout, fmt="spice")
* Anonymous `circuit.Package`
* Generated by `vlsirtools.SpiceNetlister`
*
.SUBCKT Straight__
+ e1 e2
* No parameters
rr1
+ e1 e2
+ 1000
* No parameters
.ENDS
Creating our Stimulus#
Let’s first look into how to map a numpy data into a SPICE waveform. We will create a testbench using the hdl21
interface. So the first thing is that we need to add our stimulus sources, or otherwise where would our pulses come from. We need to construct this into our testbench module. This means we need to generate the connectivity of our signal sources in relation to the ports of our circuit. We create a custom testbench module that we will use to perform this simulation. This needs to
contain also our voltage sources for whatever test that we would be performing.
[ ]:
@h.module
class OperatingPointTb:
"""# Basic Extracted Device DC Operating Point Testbench"""
VSS = h.Port() # The testbench interface: sole port VSS - GROUND
VDC = h.Vdc(dc=1)(n=VSS) # A DC voltage source
# Our component under test
dut = example_straight_resistor()
dut.e1 = VDC.p
dut.e2 = VSS
In this basic test, we can analyse the DC operating point relationship for this circuit. Now, in order to make these type of simulations more automated to run at scale, piel
provides some wrapper functions that can be parameterised with most common simulation parameters:
A Simple DC Operating Point Simulation#
[ ]:
simple_operating_point_simulation = piel.configure_operating_point_simulation(
testbench=OperatingPointTb, name="simple_operating_point_simulation"
)
simple_operating_point_simulation
Sim(tb=Module(name=OperatingPointTb), attrs=[Op(name='operating_point_tb'), Save(targ=<SaveMode.ALL: 'all'>)], name='Simulation')
We can now run the simulation using ngpsice
. Make sure you have it installed, although this will be automatic in the IIC-OSIC-TOOLS environment:
[ ]:
results = piel.run_simulation(simulation=simple_operating_point_simulation)
results
SimResult(an=[OpResult(analysis_name='operating_point_tb', data={'v(xtop.vdc_p)': 1.0, 'i(v.xtop.vvdc)': -0.001})])
We can access the data as a dictionary too:
[ ]:
results.an[0].data
{'v(xtop.vdc_p)': 1.0, 'i(v.xtop.vvdc)': -0.001}
A Simple Transient Simulation#
Let’s assume we want to simulate in time how a pulse propagates through our circuit. A resistor on its own will have a linear relationship with the pulse, which means we should see how the current changes from the input pulse in time accordingly. For clarity, we will make a new testbench, even if there are ways to combine them.
[ ]:
@h.module
class TransientTb:
"""# Basic Extracted Device DC Operating Point Testbench"""
VSS = h.Port() # The testbench interface: sole port VSS - GROUND
VPULSE = h.Vpulse(
delay=1 * h.prefix.m,
v1=-1000 * h.prefix.m,
v2=1000 * h.prefix.m,
period=100 * h.prefix.m,
rise=10 * h.prefix.m,
fall=10 * h.prefix.m,
width=75 * h.prefix.m,
)(
n=VSS
) # A configured voltage pulse source
# Our component under test
dut = example_straight_resistor()
dut.e1 = VPULSE.p
dut.e2 = VSS
Again we use a simple piel
wrapper:
[ ]:
simple_transient_simulation = piel.configure_transient_simulation(
testbench=TransientTb,
stop_time_s=200e-3,
step_time_s=1e-4,
name="simple_transient_simulation",
)
simple_transient_simulation
Sim(tb=Module(name=TransientTb), attrs=[Tran(tstop=0.1*UNIT, tstep=0.0001*UNIT, name='transient_tb')], name='Simulation')
[ ]:
piel.run_simulation(simple_transient_simulation, to_csv=True)
When you run the simulation using the piel
run_simulation
command, there is a to_csv
flag that allows us to save the data and access it afterwards. We access the transient simulation in Pandas accordingly:
[ ]:
transient_simulation_results = pd.read_csv("TransientTb.csv")
transient_simulation_results.iloc[20:40]
Unnamed: 0 |
time |
v(xtop.vpulse_p) |
s |
|
---|---|---|---|---|
20 |
20 |
0.00107 |
-0.986 |
0.000986 |
21 |
21 |
0.00115 |
-0.97 |
0.00097 |
22 |
22 |
0.00125 |
-0.95 |
0.00095 |
23 |
23 |
0.00135 |
-0.93 |
0.00093 |
24 |
24 |
0.00145 |
-0.91 |
0.00091 |
25 |
25 |
0.00155 |
-0.89 |
0.00089 |
26 |
26 |
0.00165 |
-0.87 |
0.00087 |
27 |
27 |
0.00175 |
-0.85 |
0.00085 |
28 |
28 |
0.00185 |
-0.83 |
0.00083 |
29 |
29 |
0.00195 |
-0.81 |
0.00081 |
30 |
30 |
0.00205 |
-0.79 |
0.00079 |
31 |
31 |
0.00215 |
-0.77 |
0.00077 |
32 |
32 |
0.00225 |
-0.75 |
0.00075 |
33 |
33 |
0.00235 |
-0.73 |
0.00073 |
34 |
34 |
0.00245 |
-0.71 |
0.00071 |
35 |
35 |
0.00255 |
-0.69 |
0.00069 |
36 |
36 |
0.00265 |
-0.67 |
0.00067 |
37 |
37 |
0.00275 |
-0.65 |
0.00065 |
38 |
38 |
0.00285 |
-0.63 |
0.00063 |
39 |
39 |
0.00295 |
-0.61 |
0.00061 |
We can plot our simulation data accordingly:
[ ]:
simple_transient_plot = piel.visual.plot_simple_multi_row(
data=transient_simulation_results,
x_axis_column_name="time",
row_list=[
"v(xtop.vpulse_p)",
"i(v.xtop.vvpulse)",
],
y_axis_title_list=["v(v.xtop.vvpulse)", "i(v.xtop.vvpulse)", "o4 Phase"],
)
simple_transient_plot.savefig(
"../../_static/img/examples/04_spice_cosimulation/simple_transient_plot.PNG"
)
Extracting Instantaneous Power and Resistance#
Now, we have the instantaneous current and voltage at a point in time. You can read some of the documentation in the piel
website (TODO LINK) in order to realise the fundamental relationships of this signals. In this circuit, we know that the node v(v.xtop.vvpulse)
is at the top of the resistor, and the circuit ground vss
is at the bottom. We can read the SPICE netlist to determine this (TODO add automatic schematic mapper). The current i(v.xtop.vvpulse)
flows across this
resistor according to standard Ohm’s law. So we can extract both the instantaneous power consumption and resistance just from these graphs. Note that this means that these values are subject to insignificant variations from the solver accuracy, so you might want to round to the nearest resistance values for example in some cases.
For this simple example, it is trivial, and we will demonstrate it, but for more complex circuit examples, having a generic framework is important:
[ ]:
transient_simulation_results["power(xtop.vpulse)"] = (
transient_simulation_results["v(xtop.vpulse_p)"]
* transient_simulation_results["i(v.xtop.vvpulse)"]
)
transient_simulation_results["resistance(xtop.vpulse)"] = np.round(
transient_simulation_results["v(xtop.vpulse_p)"]
/ transient_simulation_results["i(v.xtop.vvpulse)"]
)
transient_simulation_results
Note that because the resistance is constant, the power consumption should only vary when there is a change in the signal input. The resistance will always remain constant as the signal voltage and current change, as it is a physical material, and it is easily verified.
[ ]:
transient_simulation_results.iloc[20:40]
Unnamed: 0 |
time |
v(xtop.vpulse_p) |
i(v.xtop.vvpulse) |
power(xtop.vpulse) |
resistance(xtop.vpulse) |
|
---|---|---|---|---|---|---|
20 |
20 |
0.00107 |
-0.986 |
0.000986 |
-0.000972196 |
1000 |
21 |
21 |
0.00115 |
-0.97 |
0.00097 |
-0.0009409 |
1000 |
22 |
22 |
0.00125 |
-0.95 |
0.00095 |
-0.0009025 |
1000 |
23 |
23 |
0.00135 |
-0.93 |
0.00093 |
-0.0008649 |
1000 |
24 |
24 |
0.00145 |
-0.91 |
0.00091 |
-0.0008281 |
1000 |
25 |
25 |
0.00155 |
-0.89 |
0.00089 |
-0.0007921 |
1000 |
26 |
26 |
0.00165 |
-0.87 |
0.00087 |
-0.0007569 |
1000 |
27 |
27 |
0.00175 |
-0.85 |
0.00085 |
-0.0007225 |
1000 |
28 |
28 |
0.00185 |
-0.83 |
0.00083 |
-0.0006889 |
1000 |
29 |
29 |
0.00195 |
-0.81 |
0.00081 |
-0.0006561 |
1000 |
30 |
30 |
0.00205 |
-0.79 |
0.00079 |
-0.0006241 |
1000 |
31 |
31 |
0.00215 |
-0.77 |
0.00077 |
-0.0005929 |
1000 |
32 |
32 |
0.00225 |
-0.75 |
0.00075 |
-0.0005625 |
1000 |
33 |
33 |
0.00235 |
-0.73 |
0.00073 |
-0.0005329 |
1000 |
34 |
34 |
0.00245 |
-0.71 |
0.00071 |
-0.0005041 |
1000 |
35 |
35 |
0.00255 |
-0.69 |
0.00069 |
-0.0004761 |
1000 |
36 |
36 |
0.00265 |
-0.67 |
0.00067 |
-0.0004489 |
1000 |
37 |
37 |
0.00275 |
-0.65 |
0.00065 |
-0.0004225 |
1000 |
38 |
38 |
0.00285 |
-0.63 |
0.00063 |
-0.0003969 |
1000 |
39 |
39 |
0.00295 |
-0.61 |
0.00061 |
-0.0003721 |
1000 |
[ ]:
simple_transient_plot_power_resistance = piel.visual.plot_simple_multi_row(
data=transient_simulation_results,
x_axis_column_name="time",
row_list=[
"resistance(xtop.vpulse)",
"power(xtop.vpulse)",
],
y_axis_title_list=[r"resistance ($\Omega$)", r"power ($W$)"],
)
simple_transient_plot_power_resistance.savefig(
"../../_static/img/examples/04_spice_cosimulation/simple_transient_plot_power_resistance.PNG"
)
So, we have extracted the power consumption throughout time for pulses we have configured. For the whole period of the simulation, we can extract the energy consumed as the integral of power for a time differential. Note that we expect the energy consumption of this particular circuit, where the resistor is constantly drawing current, to be constantly increasing in time. We can perform a cumulative sum over our power(xtop.vpulse)
dataframe which is a discrete integral, and we can multiply
that term with the time value to determine the total energy consumption at a point in time.
[ ]:
transient_simulation_results["energy_consumed(xtop.vpulse)"] = (
transient_simulation_results["power(xtop.vpulse)"]
* transient_simulation_results["time"].diff()
).cumsum()
transient_simulation_results
[ ]:
simple_energy_consumed_plot = transient_simulation_results.plot(
x="time", y="energy_consumed(xtop.vpulse)"
)
simple_energy_consumed_plot.get_figure().savefig(
"../../_static/img/examples/04_spice_cosimulation/simple_energy_consumed_plot.PNG"
)
A full visualisation of the signal is including the cumulative energy use:
[ ]:
simple_transient_plot_full = piel.visual.plot_simple_multi_row(
data=transient_simulation_results,
x_axis_column_name="time",
row_list=[
"v(xtop.vpulse_p)",
"i(v.xtop.vvpulse)",
"resistance(xtop.vpulse)",
"power(xtop.vpulse)",
"energy_consumed(xtop.vpulse)",
],
y_axis_title_list=[r"$V$", r"$A$", r"$\Omega$", r"$W$", r"$J$"],
)
simple_transient_plot_full.savefig(
"../../_static/img/examples/04_spice_cosimulation/simple_transient_plot_full.PNG"
)
Driving our Phase Shifter#
We have demonstrated how we can extract a simple model of a resistor and create different types of SPICE
simulations. Now, let’s consider how would this affect the photonic performance in a phase-shifter context. One important aspect is that we need to create a mapping between our analogue voltage and a phase. Ideally, this should be a functional mapping.
Let’s create an arbitrary phase mapping with a function that bounds the voltage linearly. What we can do is create a function with a particular phase-power slope. This will vary depending on the thermo-optic modulator design so we will choose an arbitrary value.
[ ]:
def linear_phase_mapping_relationship(
phase_power_slope: float,
zero_power_phase: float,
):
"""
This function returns a function that maps the power applied to a particular heater resistor linearly. For
example, we might start with a minimum phase mapping of (0,0) where the units are in (Watts, Phase). If we have a ridiculous arbitrary phase_power_slope of 1rad/1W, then the points in our linear mapping would be (0,0), (1,1), (2,2), (3,3), etc. This is implemented as a lambda function that takes in a power and returns a phase. The units of the power and phase are determined by the phase_power_slope and zero_power_phase. The zero_power_phase is the phase at zero power. The phase_power_slope is the slope of the linear mapping. The units of the phase_power_slope are radians/Watt. The units of the zero_power_phase are radians. The units of the power are Watts. The units of the phase are radians.
Args:
phase_power_slope (float): The slope of the linear mapping. The units of the phase_power_slope are radians/Watt.
zero_power_phase (float): The phase at zero power. The units of the zero_power_phase are radians.
Returns:
linear_phase_mapping (function): A function that maps the power applied to a particular heater resistor linearly. The units of the power and phase are determined by the phase_power_slope and zero_power_phase. The zero_power_phase is the phase at zero power. The phase_power_slope is the slope of the linear mapping. The units of the phase_power_slope are radians/Watt. The units of the zero_power_phase are radians. The units of the power are Watts. The units of the phase are radians.
"""
def linear_phase_mapping(power_w: float) -> float:
"""
We create a linear interpolation based on the phase_power_slope. This function returns phase in radians.
"""
return phase_power_slope * power_w + zero_power_phase
return linear_phase_mapping
[ ]:
our_phase_power_map = (
piel.models.physical.electro_optic.linear_phase_mapping_relationship(
phase_power_slope=10, zero_power_phase=1
)
)
power_w_i = np.linspace(0, 1)
linear_phase_power_mapping = piel.visual.plot_simple(
x_data=power_w_i,
y_data=our_phase_power_map(power_w=power_w_i),
ylabel=r"Phase ($\phi$)",
xlabel=r"Power ($W$)",
)
linear_phase_power_mapping.savefig(
"../../_static/img/examples/04_spice_cosimulation/linear_phase_power_mapping.png"
)
This is all very nice and good, but now we need to map our power in time, to the corresponding phase of our phase shifter. We follow the same principles as the previous digitally-driven modulator example.
[ ]:
our_resistive_mzi_2x2_2x2_phase_shifter_optical_netlist = (
our_resistive_mzi_2x2_2x2_phase_shifter.get_netlist(exclude_port_types="electrical")
)
[ ]:
mzi2x2_model, mzi2x2_model_info = sax.circuit(
netlist=our_resistive_mzi_2x2_2x2_phase_shifter_optical_netlist,
models=piel.models.frequency.get_default_models(),
)
[ ]:
mzi2x2_analogue_active_unitary_array = list()
for power_i in transient_simulation_results["power(xtop.vpulse)"]:
phase_i = our_phase_power_map(power_w=power_i)
mzi2x2_active_unitary_i = piel.sax_to_s_parameters_standard_matrix(
mzi2x2_model(sxt={"active_phase_rad": phase_i}),
input_ports_order=(
"o2",
"o1",
),
)
mzi2x2_analogue_active_unitary_array.append(mzi2x2_active_unitary_i)
transient_simulation_results["unitary"] = mzi2x2_analogue_active_unitary_array
Note that this operation is very computationally intensive, as we are getting the unitary at every particular point in time within the analog simulation resolution. However, it is reasonable computable within about a minute. We can now see how our optical ports change following the principles in the digital simulation.
[ ]:
optical_port_input = np.array([1, 0])
output_amplitude_array_0 = np.array([])
output_amplitude_array_1 = np.array([])
for unitary_i in transient_simulation_results.unitary:
output_amplitude_i = np.dot(unitary_i[0], optical_port_input)
output_amplitude_array_0 = np.append(
output_amplitude_array_0, output_amplitude_i[0]
)
output_amplitude_array_1 = np.append(
output_amplitude_array_1, output_amplitude_i[1]
)
transient_simulation_results["output_amplitude_array_0"] = output_amplitude_array_0
transient_simulation_results["output_amplitude_array_1"] = output_amplitude_array_1
transient_simulation_results["output_amplitude_array_0_abs"] = np.abs(
transient_simulation_results.output_amplitude_array_0
)
transient_simulation_results["output_amplitude_array_0_phase_rad"] = np.angle(
transient_simulation_results.output_amplitude_array_0
)
transient_simulation_results["output_amplitude_array_0_phase_deg"] = np.angle(
transient_simulation_results.output_amplitude_array_0, deg=True
)
transient_simulation_results["output_amplitude_array_1_abs"] = np.abs(
transient_simulation_results.output_amplitude_array_1
)
transient_simulation_results["output_amplitude_array_1_phase_rad"] = np.angle(
transient_simulation_results.output_amplitude_array_1
)
transient_simulation_results["output_amplitude_array_1_phase_deg"] = np.angle(
transient_simulation_results.output_amplitude_array_1, deg=True
)
transient_simulation_results
Let’s plot this now:
[ ]:
simple_ideal_o3_mzi_2x2_plots = piel.visual.plot_simple_multi_row(
data=transient_simulation_results,
x_axis_column_name="time",
row_list=[
"power(xtop.vpulse)",
"output_amplitude_array_0_abs",
"output_amplitude_array_0_phase_deg",
],
y_axis_title_list=["e1 Power", "o3 Amplitude", "o3 Phase"],
)
simple_ideal_o3_mzi_2x2_plots.savefig(
"../../_static/img/examples/04_spice_cosimulation/simple_ideal_o3_mzi_2x2_plots.PNG"
)
[ ]:
simple_ideal_o4_mzi_2x2_plots = piel.visual.plot_simple_multi_row(
data=transient_simulation_results,
x_axis_column_name="time",
row_list=[
"power(xtop.vpulse)",
"output_amplitude_array_1_abs",
"output_amplitude_array_1_phase_deg",
],
y_axis_title_list=["e1 Phase", "o4 Amplitude", "o4 Phase"],
)
simple_ideal_o4_mzi_2x2_plots.savefig(
"../../_static/img/examples/04_spice_cosimulation/simple_ideal_o4_mzi_2x2_plots.PNG"
)
Model Capabilities#
We are able to arbitrarily create a model of our devices, and simulate their performance when driven by a particular circuit or not. This is useful as long as our analog model accurately represents the analog performance of our device. It is important, also, to consider opto-physical timing effects not accounted for by transient analog electronic simulations such as thermo-optic relaxation times or carrier-mobility effects. In this case, look into our component-codesign basics example.