Large-Scale Crosstalk-Corrected Thermo-Optic Phase Shifter Arrays in Silicon Photonics

We introduce a thermo-optic phase shifter (TOPS) array architecture with independent phase control of each phase shifter for large-scale and high-density photonic integrated circuits with two different control schemes: pulse amplitude modulation (PAM) and pulse width modulation (PWM). We realize a compact spiral TOPS and a 288-element high-density row-column TOPS array with this architecture and drive TOPS with waveforms of both control schemes and of different array sizes. We present a thermal excitation model and a finite difference method-based simulation to simulate large-scale TOPS arrays and compare both schemes experimentally and theoretically. We also analyze the effects of thermal crosstalk in the realized TOPS array and implement a thermal crosstalk correction algorithm with the developed model. The high-density TOPS array architecture and the thermal crosstalk correction algorithm pave the way for high-density TOPS arrays with independent phase control in large-scale photonic integrated circuits interfaced with electronics limited in voltage swing and bandwidth.


I. INTRODUCTION
I NTEGRATED photonics allows numerous bulk optic components such as lenses, modulators, and fiber-optics to be integrated on a thin substrate to realize millimeter-scale, energy-efficient optical systems with better manufacturability and cost [1]. The advent of large-scale photonic integrated circuits (PICs) largely driven by silicon photonics enabled high-performing systems driving practical applications in sensing, communications, and computing [2].
Silicon photonics leverage the high refractive index of silicon to enable the realization of optical waveguides with tight mode confinement. This allows compact nanophotonic structures to be implemented on chip at a high volume. Combined with its CMOS compatibility, silicon photonics enable the largest scale PICs with tens of thousands of components realized on a single chip [3]. Consequently, many large-scale systems such as optical phased arrays [4]- [6], programmable photonic circuits [7], optical neural networks [8], and quantum photonic processors [9] have been implemented on silicon photonics platform.
Large-scale PICs ubiquitously utilize phase shifters for modulation, tuning, calibration of systematic phase errors, and correction of random phases along the signal path [10]. As the circuits scale, a larger number of phase shifters are needed for precise signal manipulation. Thermo-optic phase shifters (TOPS) are a good candidate for high-density silicon PICs because of their low optical loss and small form factor [11], [12]. However, high-density TOPS arrays suffer from thermal crosstalk [13].
In high-density integration of TOPS, unless photonics and electronics are both monolithically integrated, electrical connections required to drive each phase shifter in a direct addressing scheme can pose a problem for scaling. Furthermore, if each phase shifter requires its own driver, the complexity of the control electronics and their power consumption scale rapidly [14], becoming a limiting factor for scalability.
To overcome these challenges, row-column or matrix addressing schemes have been implemented [4]- [6], [15]. In these schemes, the number of electrical connections scales as N with N 2 phase shifters. This reduces the complexity of the electronic interface as opposed to a direct addressing scheme where each phase shifter has a separate connection and where the number of electrical connections scales as N . In the rowcolumn architecture, each phase shifter is connected between a column and row node in series with a diode. The diode serves as a switch for the current to flow between a column and row node and modulate the phase shifter. In such a scheme, columns determine which diodes will be in forward bias, while rows determine the amplitude of the electrical signal. Electrical signals used to modulate the TOPS are time multiplexed across multiple columns, by which each TOPS receives an average electrical power proportional to the duty cycle as shown in Fig. 1. The average electrical power received by each TOPS can then be controlled by either pulse-amplitude modulation (PAM) or pulse-width modulation (PWM) of the row voltages.
In this paper, we realize a compact TOPS design and a largescale 288-element row-column TOPS array with independent phase control using a folded row-column architecture [4] along with integrated Mach-Zehnder interferometers (MZIs) and row-column photodiodes (PDs) to probe the TOPS phase shifts and calibrate the array. Driving a TOPS+MZI test structure with PAM and PWM waveforms of different array sizes, we compare PAM and PWM driving schemes. We introduce a thermal excitation model accompanied with a finite difference method-based simulation to predict the experimental results and analyze both schemes. We also showcase the effects of thermal crosstalk between phase shifters and implement a thermal crosstalk correction algorithm to completely cancel the crosstalk in the row-column TOPS array.
where λ 0 is the wavelength in free space and L is the length of the heated waveguide. Then, the phase shift in the TOPS is where ∆T is the temperature difference, and γ is the temperature coefficient of refractive index. Then, ∆T π , temperature difference for π phase shift is A single TOPS can be modeled as a simple RC circuit, with P π , power required for π phase shift, defined as where A is the effective surface area and k is the effective thermal conductivity between the waveguide and the heat sink. In the RC model, G = 1/R = Ak is the thermal conductance of TOPS to the heat sink. Another parameter to characterize a TOPS is its time constant, τ , similarly defined with the first-order RC approximation as τ = C Ak (5) where C is the heat capacity of TOPS.

A. Thermal Excitation Model
We now develop a time-varying lumped model starting with a 3D analysis of TOPS seen in Fig. 2, which shows the cross-section of the realized TOPS design and array. Without loss of generality, we consider a heater on one side of the waveguide in an XYZ space. Here, we analyze the dynamics between four structures: heater, waveguide, substrate, and the environment acting as a perfect heat sink. We assume uniform material for all structures, uniform temperature distribution within each structure and in X and Z directions. We define a time-varying temperature function, T (y, t), that characterizes the local instantaneous temperature. When some instantaneous electrical power, P (t), is applied to the resistive heater, this power is expended through four mechanisms (ignoring power dissipated in forms other than Joule heating): • Heat used to increase the heater temperature (q H = ρ H c H dT (y, t)dV ) where ρ H , c H and V are the heater mass density, specific heat capacity and volume.) • Heat exchange with the waveguide (q W = −k W H dtd A · ∇T (y, t) where k W H is the effective thermal conductivity between the waveguide and the heater, d A is the differential surface area vector, and ∇T (y, t) is the instantaneous temperature gradient.) • Heat exchange with the substrate (q S = −k SH dtd A · ∇T (y, t) where k SH is the effective thermal conductivity between the substrate and the heater.) • Heat exchange with surrounding environment other than the waveguide and substrate (q E = −k EH dtd A·∇T (y, t) where k EH is the effective thermal conductivity between the environment and the heater.) 1) No Crosstalk: First, we neglect the crosstalk between phase shifters. Then, dP (t)dt = q H + q W + q S + q E . Going through the derivation in the Supplementary Material, we have the matrix ODE for the lumped model of the system.
where we defined G ij = 1/R ij = A ij k ij . Since we assumed uniform temperature distribution within each structure, k ij reduces to effective thermal contact conductivity between the respective structures. Furthermore, D H = G W H +G SH +G EH , D W = G W H +G SW +G EW , and D S = G SH +G SW +G ES are the sum of thermal conductances for the matrix diagonal elements. Note that T E is a constant as the environment acts as a heat sink. We also take note that substrate usually facilitates the thermal conduction between structures. Therefore, 2) With Crosstalk: Now, we include the crosstalk between phase shifters and repeat the analysis (see Supplementary Material) by scaling the single TOPS model to N phase shifters. Hence, we have the generalized matrix ODE.
We also again take note that substrate facilitates the thermal conduction between structures. Therefore, for each TOPS, This also means crosstalk between TOPS is entirely due to G S,iS,j = 1/R S,iS,j . Hence, 1 We define D H,i = S H,iH + S W H,i + S SH,i + G EH,i , For any array of phase shifters, this system of ODEs can be solved numerically or through variation of parameters to obtain the time-domain response of the heaters, waveguides, and substrate. This enables fast modeling of large-scale TOPS arrays without the need for computationally-intensive FDTD simulations enabling efficient engineering of driving waveforms for arbitrary system parameters.

B. Finite Difference Simulation
Using the framework developed with the thermal excitation model, we implement a finite difference method-based timedomain simulation, numerically solving the aforementioned matrix ODEs to simulate large-scale TOPS arrays. Electrical signals to drive the TOPS array are defined by P i . Therefore, we construct an ideal signal and pass it through an FIR low pass filter designed to simulate the bandwidth of the electronics. Then, to reduce the matrix size and decrease computation time in the simulation, we make some simplifications to the model. We make a first-order approximation for the frequency range we are characterizing the array. We set the substrate parameters, C S and G ES significantly high, making the substrate an ideal heat sink. Since substrate heating introduces a low-frequency pole to the system, this simplification makes the simulation at low frequencies less accurate but can be removed from the simulation. Furthermore, we define the heater to have low C H , which is a valid assumption for our design as will be shown by the first-order frequency response in our characterization. However, for different designs, this simplification can also be modified. For thermal crosstalk simulations, we assume only first-order crosstalk, meaning G S,iS,i±n = G (n) SS = 0 with n > 1. We finally approximate G = G SW = G SH . With these simplifications, we rederive (3), (4), and (5) to find the remaining simulation parameters, namely C W , G, G SS , and MZI extinction ratio (extinction ratio due to non-ideal couplers), ER (see Supplementary Material).
where P (1) π is the power required to drive a TOPS adjacent to the main TOPS to induce π phase shift at the main TOPS and P min , P max correspond to the observed minimum and maximum optical powers, respectively. C W , G, G SS , and ER can all be characterized experimentally.

III. ROW-COLUMN TOPS ARRAYS
A row-column TOPS array can be realized by arranging multiple TOPS with diode switches in a row-column format. It can be driven with PAM [5], [6] or PWM [4], [15], [16] as seen in Fig. 3. The average AC power dissipated (P AC ) is determined by the row voltage amplitude for PAM and by the row duty cycle for PWM. DC power can be calculated from row and column voltages via P DC = ∆V 2 /R, where R is the heater resistance. Also note that for PWM, there are additional terms (i.e. P Φ , V Φ ) that correspond to the minimum pulse width (t Φ ) resolvable by the limited bandwidth of the row driver.
We realize a 32 × 9 row-column TOPS array on a silicon photonics platform as seen in Fig. 4. The row-column array utilizes a folded architecture, fitting two electrical rows per one physical row with a shorter optical path length, reducing the optical loss and increasing the coherence between branches. This architecture also reduces the electrical length of columns by half that allows a smaller metal width for routing, improving scalability.
Light is coupled to the chip via a grating coupler with 3.5 dB insertion loss and is distributed to the TOPS array via a 1:256 splitter tree with 2.4 dB excess loss. There are 32 heaters without any optical waveguides acting as dummy phase shifters [4]. Each TOPS is comprised of a 300 µm spiral waveguide surrounded by doped silicon serving as a Joule heater with an area of 1320 µm 2 . The heater is fitted with a silicon diode in series to switch on/off the TOPS in the array. There is also an oxide etch around each TOPS to minimize thermal crosstalk. TOPS is designed to have, on average, a P π of 10 mW, a bandwidth of 30 kHz, and an insertion loss of <0.2 dB. To characterize the phase shifts from TOPS and use them as feedback information, we also realize an MZI with each TOPS pair that taps the TOPS outputs with 95:5 couplers and combines them with 0.3 dB excess loss. The output is fed to a PD in a row-column PD array. The number of electrical connections for the PD array is also reduced from N 2 to N with the same row-column addressing scheme by adding a silicon diode switch in series with the PD [17]. The total die area with the splitter tree, MZI+PD array, and electrical pads is 6 mm 2 . The total excess loss for each channel is 6.5 dB. Circuit diagrams of the TOPS array are seen in Fig. 5.
We characterize a TOPS test structure for the simulation parameters as shown in Fig. 6. We find its P π to be 10.6 mW and its bandwidth to be 56.8 kHz. These parameters correspond to C W = 440 pJ/K and G = G SW = G SH = 62.8 µW/K, which are within an order of magnitude of the expected parameters from material properties given 0.5 µm waveguide width, 300 µm waveguide length, 220 nm waveguide height, and 2 µm BOX layer height (C W = 54.7 pJ/K with Si specific heat capacity of 20 J/mol · K [18] and G = 71 µW/K with SiO 2 thermal conductivity of 1.1 W/mK [19]). To characterize crosstalk, we drive a dummy phase shifter adjacent to the main TOPS without driving the main TOPS and find a P Due to fabrication variations, the characterized values from this parameter extraction can slightly differ from one TOPS to another in the array. Furthermore, TOPS placement in the array and thermal properties of the experimental setup may change these numbers. Therefore, we usually use the characterized parameters or simulated design specifications only as a starting point to fit the simulation parameters to the measured data. We also note that, as seen in Fig. 6b, the single TOPS frequency  response is of first order, making the first-order approximation in our simulation valid within the frequency range we take measurements.
For an arbitrary array with M rows and N columns such as in Fig. 3, there are M × N phase shifters. For PAM, the duty cycle of each phase shifter is fixed at 1/N . For PWM, the duty cycle of each phase shifter is varying but has a maximum of 1/N . We also define a DC power (P DC ) that corresponds to the average AC power dissipated in the TOPS (P AC ). Then, for PAM, the DC power is varying but has a maximum of N × P 2π , where P 2π is the phase shifter power required for 2π phase shift. For PWM, the DC power is fixed at N × P 2π . Therefore, P DC = P AC /D for both PAM and PWM, where D is the duty cycle.

A. Control Electronics
We first drive a single TOPS test structure with different waveforms to simulate row-column TOPS arrays of different sizes without thermal crosstalk and then the 32 × 9 TOPS array to investigate and correct thermal crosstalk. The control electronics for these measurements are shown in Fig. 7. The PIC with the TOPS array is wirebonded to an interposer board sitting on top of a custom-made motherboard with bias circuitry and TIAs for 16 PD channels (16×16 row-column multiplexed channels on-chip) and routing circuitry for driving 32 rows and 9 columns. PD channel multiplexing is achieved via a 4-to-16 demultiplexer, whose outputs are used to bias the PD rows, whereas PD columns are connected to TIAs for read-out. TIA outputs are digitized with 16 channels of 100 MSa/s, 14-bit digitizers. The switching columns are driven by a custom-made board with a 4-to-16 demultiplexer and operational amplifiers for 9 channels. The driving rows are driven by 32 channels of 1 GSa/s, 14-bit arbitrary waveform generators (AWGs) with an on-board FPGA. Row and column drivers share the same ground and are time-synchronized. Time synchronization is accomplished by a 4-bit synchronous counter with reset functionality that drives the logic inputs for the column driving demultiplexer. Desired waveforms for the 32 rows are pre-loaded to the row driver AWGs, and a separate waveform acting as the clock for the counter is sent to the counter board. Between each cycle, AWG also sends a trigger signal to the counter's reset pin to reset the switching cycle. To reduce noise coupled from the digital circuitry, driving and read-out electronics have two separate grounds with appropriate filtering for all the DC power supply.
We now investigate both PAM and PWM driving schemes in more detail.

B. Pulse Amplitude Modulation
PAM driving entails controlling the average power dissipated in the heaters by tuning the amplitude of the square waveform that has a constant pulse width as observed by a single TOPS. First, without crosstalk, we drive a single TOPS (P π = 13.4 mW, τ = 20.9 µs, ER = −21.3 dB) with this waveform at different frequencies for N = 2 at P π as seen in Fig 8. The electronics bandwidth of the electronics is 200 MHz, which determines the cut-off frequency of the low-pass filter in the simulation.
We then drive another single TOPS (P π = 11.6 mW, τ = 20 µs, ER = −19.9 dB) at different frequencies with different AC powers to compare the simulated and experimental MZI curves seen in Fig. 9. We note that at frequencies above 1 MHz, extinction ratio measurement of TOPS is dominated by ER, extinction ratio same as that of a directly-addressed array. Since TOPS is driven at 2P π DC power with 1/2 duty cycle, there is an observed phase swing around π. This phase swing is reduced at higher frequencies, leading to a higher extinction ratio. The discrepancy between simulation and measurement at 1 kHz is due to the first approximation, through which we set C S and G ES significantly high, making the substrate a perfect heat sink. However, in measurement, substrate temperature rises a measurable amount at high AC powers.

C. Pulse Width Modulation
PWM driving entails controlling the average power dissipated in the heaters by tuning the duty cycle of the square waveform that has a constant amplitude as observed by a single TOPS. Without crosstalk, we again drive a single TOPS (P π = 14.1 mW, τ = 23.1 µs, ER = −21.3 dB) with this waveform at different frequencies for N = 2 at P π as seen in Fig. 10.
We then again drive another single TOPS (P π = 14.3 mW, τ = 18.6 µs, ER = −23.1 dB) at different frequencies with different AC powers to compare simulated and experimental MZI curves as seen in Fig. 11.

D. PAM vs. PWM
With these measurements, we can compare PAM and PWM, whose waveforms and some of the definitions that will be For PWM, since TOPS is driven at 4P π DC power with 1/4 duty cycle, there is a phase swing around π with a greater reduction in extinction ratio at lower driving frequencies compared to PAM.
used in the following analysis are shown in Fig. 3. PWM driving offers benefits in linearizing the control scheme and simplifying the electronics since voltage does not need to be tunable, whereas PAM requires tunability in amplitude. However, leveraging the thermal excitation model, a more detailed analysis shows more striking differences between PAM and PWM with respect to driving electronics specifications. Hence, we compare the two schemes regarding electronics' voltage swing and bandwidth.

1) Voltage Swing:
The voltage swing requirement for either scheme is set by the required phase shift range, which in most cases is 2π for full-wave coverage. For a single TOPS (N =1), the voltage swing requirement for PAM is V 2π , and P 2π = V 2 2π /R, where R is the TOPS heater's electrical resistance. However, for PWM, since AC power is only tunable with duty cycle, an AC power that requires a pulse width corresponding to a bandwidth significantly larger than the electronics' bandwidth cannot be realized accurately. This can be observed in Fig. 11 at low AC powers (<2 mW). To realize phase shifts corresponding to these AC powers, voltage swing requirement needs to be V 2 2π + V 2 Φ , where V Φ corresponds to the AC power required to recover the lost small phase shifts (P Φ ). This couples the bandwidth and voltage swing requirements for PWM.
Moreover, the voltage swing requirement scales linearly with √ N , and for full-wave coverage, the driving electronics' output voltage swing should be √ N V 2π for PAM and N (V 2 2π + V 2 Φ ) for PWM. Therefore, the PAM voltage swing requirement is less stringent than that of PWM especially for large-scale TOPS arrays if P Φ > P DAC where P Φ is the power corresponding to t Φ and P DAC is the power corresponding to t DAC , the sampling resolution of the row driver. This increased voltage swing requirement may become significantly detrimental to large-scale TOPS arrays with PWM driving since it also implies a larger reverse breakdown voltage requirement for the TOPS diodes and could extend the peak voltage requirement beyond the CMOS transistor breakdown voltage [4], [14]. 2) Bandwidth: It can already be seen from the amplitudedomain responses that PAM driving reaches a higher extinction ratio at lower frequencies than PWM. To demonstrate this more clearly, we do a frequency-domain response comparison by finding the extinction ratio of the amplitude-domain response at each frequency for N = 2. The resulting response obtained both experimentally and with the simulation is shown in Fig. 12. As seen here, for systems where bandwidth is the primary concern, PAM driving offers the benefit of reaching a higher extinction ratio than PWM driving. The reason for this can be explained by investigating Figs. 8 and 10 more carefully.
In the time domain, TOPS phase shift oscillates back and forth centered around the intended phase shift. Phase swing of this oscillation determines the phase precision/extinction ratio of the TOPS/MZI. Phase swing is positively correlated with the time when TOPS is not driven, t OFF , the steadystate temperature difference of the waveguide when TOPS is driven vs. not driven, ∆T , and is negatively correlated with the time constant, τ , as determined by the solution to the system of ODEs in the model. While τ is specific to the TOPS design and ∆T is dictated by the required phase shift, t OFF increases with pulse width, which is determined by electronics' bandwidth. PAM has an inherent advantage since its t OFF is lower for a certain phase shift than that of PWM, t OFF,PAM ≤ t OFF,PWM . The difference between these times, ∆t OFF = t OFF,PWM − t OFF,PAM , determines the difference in the extinction ratio of these two schemes. Therefore, as the modulation frequency increases, ∆t OFF decreases, consequently reducing the difference in extinction ratio. The same effect takes place as N is scaled, meaning both PAM and PWM bandwidth requirement scales linearly with N . To demonstrate this, we plot the difference in extinction ratio by performing the frequency-domain analysis in Fig. 13 for different N . As seen in Fig. 13, difference in the extinction ratio between PAM and PWM is minimal for large N and for high bandwidth. However, this difference is significant (up to 4.7 dB) for bandwidth-limited systems with small N TOPS arrays.
Positive correlation between phase swing and ∆T also explains the reduction in the extinction ratio as power dissipation increases for higher phase shifts (e.g. 2π) in Figs. 9 and 11.
Finally, since voltage swing and bandwidth requirements are coupled for PWM, bandwidth could be used to minimize the voltage swing requirement by setting P Φ < P DAC , namely t Φ < t DAC . However, since t Φ ∝ 1 BW , as we minimize the voltage swing requirement for PWM, bandwidth requirement becomes more stringent. Therefore, PWM suffers from a more stringent voltage swing requirement and/or a more stringent bandwidth requirement.

IV. THERMAL CROSSTALK
Now, we showcase the effects of thermal crosstalk in rowcolumn TOPS arrays and demonstrate a correction algorithm with the information fed back from the integrated MZI+PD array and utilizing the previously detailed thermal excitation model. To maximize the effects of crosstalk, we select a TOPS in the middle of the array (row 16, column 5). To characterize the phase shift from the TOPS with and without crosstalk, we use the TOPS left of the main TOPS as the reference for the corresponding MZI. First, we drive the main TOPS with N = 9 PAM of 200 kHz frequency (0.56 µs pulse width) but still keep the other TOPS off to have a baseline for the phase shift without crosstalk. Then, we drive 142 TOPS with the same waveform. We only keep the reference TOPS off since it serves as a reference and the TOPS right of the main TOPS off to introduce an asymmetry to observe the crosstalk effect. We drive the rest of the TOPS (16 rows x 9 columns) with a maximum amplitude of 9 V to maximize crosstalk. We also simulate the amplitude-domain response with the thermal excitation model. We use the simulation parameters (C W =440 pJ/K, G=62.8 µW/K, G (1) SS =4.57 µW/K, ER = −0.665 dB) extracted from the TOPS characterization as a starting point. We tune G (1) SS = 4.57 µW/K to G (1) SS = 13.2 µW/K since we neglected higher-order crosstalk in the simulation and crosstalk dynamics in the array are expected to be different from the test structure and from setup to setup. Hence, measured and simulated amplitude-domain responses for the main TOPS with and without crosstalk are shown in Fig. 14 with linear optical power on the y-axis. The measurement has a low extinction ratio due to stray light in the substrate reaching the integrated PD array but can be improved by adding Ge between components that absorbs stray light. As seen here, thermal crosstalk causes a shift in the MZI curve. This shift consequently leads to a decrease in extinction ratio due to the aforementioned positive correlation between phase swing and ∆T .

A. Thermal Crosstalk Correction
We now demonstrate a thermal crosstalk correction algorithm utilizing the thermal excitation model. For TOPS arrays with direct addressing, various thermal crosstalk correction algorithms have been proposed based on eigenmode decomposition [20], [21] and matrix inversion [22]. These algorithms can also be implemented for row-column TOPS arrays to map the phase shifts from driving with crosstalk to driving without crosstalk. However, due to the phase swing caused by temporal multiplexing, the change in extinction ratio also needs to be corrected. In conjunction with our model, we implement an algorithm based on matrix inversion to correct thermal crosstalk. Assuming only first-order crosstalk and a linear system ( ∂φ ∂P = constant), we define a square matrix, C, to specify the modulation and coupling coefficients at steady state. Then, where C ii+1 = Γ G (1) given Γ = 2πLγ λ0 and C ii+n = 0 given n > 1. Then, by finding the inverse of C, required P can be calculated, P = C −1 Φ. To correct the extinction ratio difference, we use the thermal excitation model to find the correct driving frequency. Simulation predicts a corrected driving frequency of 367 kHz, setting the PAM driving frequency of the crosstalk corrected measurement. Measured and simulated MZI curves throughout the crosstalk correction algorithm are shown in Fig. 15.
Since we used the parameters from the parameter extraction, this shows how our algorithm can be used in any arbitrary row-column TOPS array using only the characterization values from its test structure with minimal tuning. This also shows that in the case where parameter extraction is inaccurate, a row-column array of integrated feedback PDs and MZIs can be used to find the corrected parameters of all orders for accurate matrix inversion and extinction ratio correction.

V. CONCLUSION
We realize a 32 × 9 TOPS array based on a folded row-column array architecture with PAM and PWM control schemes and a thermal excitation model to analyze and compare both schemes. By leveraging our model, we shine light on the design requirements for the driving electronics and observe the time-domain, amplitude-domain, and frequency-domain behavior of row-column TOPS arrays. This enables waveform engineering for more versatile row-column TOPS arrays with less stringent electronic voltage swing and bandwidth requirements while maintaining the scalability advantage (M+N drivers for M × N TOPS).
We also investigate the effects of thermal crosstalk in rowcolumn TOPS arrays and implement a crosstalk correction algorithm to completely cancel the thermal crosstalk for a TOPS in the middle of our array. This demonstrates negating the effects of thermal crosstalk in row-column TOPS architectures, unlocking higher-density TOPS integration.
From PIC dimensions, it is estimated that up to more than 40,000 crosstalk-corrected TOPS can be realized on a single 30 mm × 30 mm reticle by only modularizing this design with a potential for further scaling with larger N row-column TOPS and feedback MZI+PD arrays. To the best of our knowledge, this realization signifies the largest row-column TOPS array with independent phase control ever demonstrated, paving the way for larger-scale and higher-density error-corrected PICs.   For no crosstalk, we have one TOPS with a heater, waveguide and substrate as shown in Fig. 16. With the definitions detailed in Section IIA, we first write a PDE for the heater.
Since temperature gradient is uniform over the surface area, d A · ∇T (y, t) = dA∇T (y, t). Then, Integrating over volume, where we used dV = dxdydz, and dA = dxdz. Then, where C H is the heat capacity of the heater. Simplifying, where we took G ij = A ij k ij and noted that A and k are the same physical parameters as in (5). (18) describes the process through which the applied electrical power is expended to increase the temperature of the heater, which then drives the heat exchange between the heater, the waveguide, and the substrate. The power exchanged with the waveguide is again expended through three mechanisms (again ignoring dissipated power in forms other than Joule heating): 1) Heat used to increase the waveguide temperature 2) Heat exchange with the substrate 3) Heat exchange with surrounding environment other than the substrate Then, where C W is the heat capacity of the waveguide and G SW , G EW signifies the thermal coupling between the waveguide and the substrate, the environment respectively. Similarly, we can also derive a PDE for the substrate.
where C S is the heat capacity of the substrate and G ES signifies the thermal coupling between the substrate and the environment. Since we assumed uniform temperature distribution within all structures in all directions, spatial temperature distribution within each structure can be ignored, simplifying spatial partial derivatives to a lumped model. This is due to (for a linearly-varying spatial temperature distribution within each structure) where ∆T i is the temperature difference within structure i, ∆T ij is the temperature difference between structures i and j at their interface, k i is the thermal conductivity of the structure i, and k ij is the thermal contact conductivity between the structures i and j. Considering only the heat flow at the interfaces between structures, G ij reduces to A ij k ij , where k ij is now the effective thermal contact conductivity between the respective structures. Defining T H , T W , T S , and T E as lumped time-varying temperatures for the heater, waveguide, substrate and environment, respectively, we simplify (18) to Similarly, (19) becomes and (20) becomes Hence, using using (22), (23), and (24), we write the matrix ODE in (6) for this lumped system.

B. With Crosstalk
We now redo the above derivation including the crosstalk between all phase shifters. If we have N phase shifters, for phase shifter i, (18) becomes To find the equations for G and C W , we draw the RC circuit model corresponding to the matrix ODE in (6) as shown in Fig. 17. With the approximations introduced in Section IIB (G ES G SW = G SH = G and C S C W C H ), the reduced circuit model becomes as shown in Fig. 18. Neglecting C W for DC analysis, ∆T π at the waveguide node can be calculated as Rearranging (31), Now, with C W , ∆T π = P π R eq R eq + 2R where R eq = R 2 1 jωC W = R 2+jωRC W . Then, ∆T π = P π R 5R + 2jωR 2 C W = P π 1 5 + j2ωRC W To extract G SS , we draw the circuit diagram for a coupled TOPS pair of first-order crosstalk, again with the same approximations, as shown in Fig. 19. Solving for the temperature at W (1) gives the following transcendental equation.