## Abstract

When solving, modeling or reasoning about complex problems, it is usually convenient to use the knowledge of a parallel physical system for representing it. This is the case of lumped-circuit abstraction, which can be used for representing mechanical and acoustic systems, thermal and heat-diffusion problems and in general partial differential equations. Integrated photonic platforms hold the prospective to perform signal processing and analog computing inherently, by mapping into hardware specific operations which relies on the wave-nature of their signals, without trusting on logic gates and digital states like electronics. Here, we argue that in absence of a straightforward parallelism a homomorphism can be induced. We introduce a photonic platform capable of mimicking Kirchhoff’s law in photonics and used as node of a finite difference mesh for solving partial differential equation using monochromatic light in the telecommunication wavelength. Our approach experimentally demonstrates an arbitrary set of boundary conditions, generating a one-shot discrete solution of a Laplace partial differential equation, with an accuracy above 95% with respect to commercial solvers. Our photonic engine can provide a route to achieve chip-scale, fast (10 s of ps), and integrable reprogrammable accelerators for the next generation hybrid high-performance computing.

### Summary

A photonic integrated platform which can mimic Kirchhoff’s law in photonics is used for approximately solve partial differential equations noniteratively using light, with high throughput and low-energy levels.

## 1 Introduction

Photonic integrated circuits (PICs) do not exist. Even if this statement could seem outrageously contradictory, we invite the reader to bear with us while we unravel the assertion. The concept of a circuit originates from connecting electronic components into a functional unit and as such is governed by certain physical rules. These fundamental rules of circuitry that, in fact, does not exist in optics and hence also not in photonic platforms. As such, the perception of a circuit applied to photonics is actually only a rather loose sense, with significant physical and technological consequences.

In electronics, circuits are simple loops, in which the flow of electrons circulates, ruled by the conservation laws, governed by quasi-static approximation of Maxwell’s equation, resulting in Kirchhoff’s law. In contrast, in photonics, light does not have return loops but is usually conveyed to be ultimately detected. Conservation laws still hold considering the light dissipation and transition into other domains, but the “flow of photons” follows Maxwell’s.

Another main aspect to consider, as a consequence of the quasi-static assumption, is that we usually refer to circuits when we can approximate their components as concentrated at singular points in space (“lumped circuits”) in which the physical quantities, such as potentials and currents, are function of time only. This approximation is possible because the wavelength of the signals, and their time-scale variation, is significantly longer and slower than the physical dimension and variations of the circuit itself, respectively. The consequence of this approximation is nonlocal effects, i.e., elements of the circuits are coupled, and local variations will affect the global performance of the circuit.

In photonics, that is not the case (Figure 1 i–iii); photonic platforms have to be considered as distributed networks since they are characterized by a footprint that is several orders of magnitude larger compared to their operating wavelength (Length_{PIC} >> *λ*_{NIR}). Typically, when designing “photonic circuits,” the aim is to design components that will generate photons and efficiently convey them, modulated in some prearranged way, to obtain a certain functionality. Thus, photonic circuits are information pathways rather than circuits (Figure 1 iii). For this reason, photonics is considerably more related to data transport than to algorithmic transformations and operations.

The physical origin of the difference between the behavior displayed by electronics versus optical circuits rests in the vastly different momenta of electronics in matter and photons of the near IR and visible frequency ranges.

Therefore, strictly speaking, the “photonic circuit”, which confines and manipulates light, effectively performs its functionality on a component-by-component manner. However, Kirchhoff’s law and nonlocal effects contribute in mapping in integral way fundamental mechanism used for describing important phenomena directly into electrical circuits and their algorithmic functionality. Indeed, such “coupled-nodes” as part of a phase-constant network can be a powerful tool for mimicking integral-differential problems, for instance, which ubiquitously pertains to a plethora of diverse scientific and engineering problems, onto photonic hardware, thus potentially relevant for future analog computing accelerators.

To further explain this concept, at the dawn of the digital era, analog processors based on electrical mesh and Kirchhoff’s law have been conceptualized and demonstrated [1], [2], [3]. Such analog processors were able to solve second order partial differential equations (PDEs) using finite difference methods (FDM), relying on continuous signals and programmed by changing the interactions between its computing elements, e.g., impedances, using minimum stored programs or algorithms; thus obtaining solutions in a completely asynchronous manner, providing one-time, (noniterative) computations independent, at first order, of the problem complexity. However, the complexities of an effective integration of a high-speed programmable and energy-efficient static-like analog mesh and the concurrent advancement of digital electronics architectures, eclipse this Kirchhoff’s electrical FDM approach.

Current (von Neumann) processors solve PDEs through numerical methods, involving iterative high precision vector-matrix precision operations, which can be both power and time-costly according to the complexity and resolution of the problem. These bottlenecks are only softened by parallel hardware (i.e., multicore processing), which do not offer a significantly different path to accelerate PDEs, due to the parallelism overhead and disadvantageous speed-up scaling with respect to the number of processors [4].

Since the late 2000s, the computing paradigm has shifted again; it seems apparent that a new class of hybrid hardware is emerging, i.e., co-processor and accelerator, able to perform a task efficiently, by homomorphically mapping a specific problem category on application-specific hardware, which solves the problem at hand in an entirely parallel manner. Such outlook can potentially ameliorate the computing pressure on digital electronics.

In this regards, integrated photonics-based signal processing, thanks to the electromagnetic nature of its signals and availing their efficient interactions with matter, places itself as a compelling solution for optical communication [5], quantum information processing [6], optical computing [7] and especially neuromorphic computing with remarkably reduced energy consumption and accelerate intelligence prediction tasks [8], [9], [10], [11]. Recently, inverse-designed metamaterial platform interfaced with integrated photonics showed the possibility of solving integral equations using monochromatic electromagnetic radiation [12]. In the field of PDE solvers, all-optical reconfigurable module based on microring resonators can solve ordinary first and second-order temporal differential-equation [13], [14].

Here we argue and demonstrate that, to a certain extent, we can induce circuit-related homomorphism by imposing an integrated photonics platform to behave as a circuit in which Kirchhoff’s law can be re-implemented. We experimentally demonstrate a photonic node, termed Kirchhoff’s photonic node (KPN), which splits the incoming light evenly into three remaining directions with minimized reflections, which can be reprogrammed emulating an optical version of Kirchhoff’s current law using optical light intensity. When multiple nodes are interconnected and arranged to generate a two-dimensional mesh, as envisioned by Anderson et al. and Sorger et al. [15], [16], in analogy to a uniform electronic resistive circuit, we are able to approximately solve a stationary partial differential equation which simulates in a noniterative fashion the heat-transfer problem of a film, via an optical finite difference method (FDM). We observe that when the proposed photonic circuit, here termed Silicon Photonic Approximate Computing Engine (SPACE), is forced to pseudo-homomorphically map the PDE, we can achieve solution accuracy up to 97% compared to a simulated heat transfer problem with the same mesh resolution which can be obtained in just 16 ps.

## 2 Results

### 2.1 Photonic Kirchhoff’s node

The fundamental unit of an electrical circuits are the nodes, which represent the terminals in which two or more circuit elements meet, and their distribution and interconnection determine circuit functionalities and operating conditions. At the node of electrical circuit, Kirchhoff’s current law applies; which states that current flowing into a node (or a junction) must be equal to current flowing out of it. This is a consequence of charge conservation (energy), if we can concurrently assume that the instantaneous variation of the magnetic flux outside a conductor and the change of charge in the conductor is zero (steady-state conditions). In photonics, due to the distributed nature (*l* ≥ *λ*_{eff}) of the platform used, the behavior of the electromagnetic radiation is governed by Maxwell’s equation in time-variant conditions, in which the electric and magnetic field are both function of time and position (Figure 1). Hence, for our analysis we consider as fundamental quantity of the photonic circuit the optical power instead of electromagnetic field intensity. Although, this is not sufficient for ensuring that the electromagnetic radiation flowing into a node is necessarily equal. We must point out that an equal splitting ratio across the outgoing ports in fundamentally unattainable, and not only impractical due to “nonidealities” given by (1) impedance mismatch at the discontinuity (i.e., joint between different network elements), which causes reflections, thus producing interference with the incoming radiation, (2) light scattering at the abrupt discontinuities, and (3) optical losses due to mode dispersion. The fundamental nonequal splitting ratio of the herein discussed passive, reciprocal, linear, optical four-port splitter is due to the fact that one cannot have zero reflection in all of its ports, unless at least four additional off-diagonal elements in its S-matrix are also zero. Nonetheless, prototyping a photonic node one can achieve a design where nonidealities, such as reflections and optical losses, are minimized. As a first step towards this goal, we aim to mimic an electric node (Figure 2A.i) with equal resistors onto a photonic platform, which evenly partitions the optical power flowing from one of its side.

Towards approaching an near-equal light splitting condition, this node design needs to meet the following three criteria: (a) to be symmetrical to both *x*- and *y*-axis in order to physically build the scalable optical mesh, and needs to provide a 1-to-3 equal splitting ratio; (b) the splitter needs to have tolerance to the fabrication variance since cascading the node will amplify the device variance; (c) the segment for light coupling should have the potential to be further integrated with tuning mechanisms (e.g., electro-optic means) in order to ensure reconfigurability and compute-programmability. Independently of the photonic node typology, the optical loss along each light path can be used as an equivalent resistances *R*_{i} in the electrical model.

For this aim, we used two distinctive approaches; Firstly, we follow a heuristic approach in which to obtain even splitting of the optical power using waveguide crossing assisted by four directional couplers integrated with ring resonators (Figure 2A.ii). To achieve the design of the KPN we use a heuristic process to obtain 1–3 equal power splitter in photonics by iteratively optimizing the splitting ratio using 3D full-wave numerical simulations. The resulting design comprises of four water-drop shaped rings placed close to two perpendicularly crossed waveguides to couple part of the light coming from one direction into both, the other two perpendicular directions and still let the remaining light pass through to the opposite port. Instead of using circular rings, the segments close to the straight waveguides are flattened to form a three-waveguide directional coupler. We used directional couplers to couple into the four feedback loops, and refrain from using neither perfectly circular rings and nor high-quality factor cavities to widen the spectral (and thermal) operating window such as to not having to use tuning (e.g., thermal, electro-optic) to control its resonance. In addition, a four-way waveguide crossing is the center of each Photonic Kirchhoff’s node to reduce the scattering and crosstalk at the intersection.

Secondly, an optimized inverse design approach is used by setting the design area to 5 µm and the even splitting functionality in the cardinal directions according to the following cost function to optimize

where *T* is the target transmittance and *n,* assuming the following constraints:

In the full space of fabricable devices, the optimization algorithm finds a structure (Figure 2A.ii) that meets these requirements (further details of the electromagnetic characterization of the inverse design KPN can be found in the Supplementary material, Figure S2).

The inversely designed KPN is one among the infinite number of configurations that would satisfy the optimization of the cost function. Even if characterized by a compact size (2 × 2 µm), it is significantly more intricate to fabricate with the same high yield as regular photonics due to its limited size. Additionally, this type of node would require a completely different configuration and related optimization process for mimicking different “optical power” partitioning, while the heuristic solution can be straightforwardly reconfigured by actively tuning the coupling coefficients with the ring resonators. It is worth noticing that the inverse design of a KPN is symmetric, linear, passive, and free of magnetic poles. Thus, the S-parameter matrix is equal to its transpose (

After optimizing the bending radii of the water-drop shaped rings, flattened coupler length and the gap between the ring and the straight waveguide, the splitting ratio can be tuned to 22, 23 and 22% with 12% reflection (Figure 2B) based on full-wave simulation (Lumerical 3D FDTD). Here the reflection is mainly caused by the return couplings from the three-waveguide couplers at the perpendicular ports (i.e., Output port 1 and 3). Instead of completely coupling to the perpendicular port, the light coming from the first two rings will be partially leaked to the rings on the other side and route the signal back to the input port. In terms of the fabrication process, 2% hydrogen silsesquioxane (HSQ) electron-beam resist is used due to its fine resolution and edge contrast, with isotropic dry etching process (SF_{6} and C_{4}F_{8}) to get the uniform height profile and vertical sidewall profile. More details related to the fabrication are given in Section 4.

In order to test the fabricated device, a 1550 nm continuous wave (CW) laser is used as light source which is coupled to the waveguide by means of a periodical grating coupler (details regarding the grating coupler are discussed in the Supplementary material). To read out the output values at each port of the KPN in parallel, an InGaAs infrared camera (Xenics) is used to capture the microscope image of the light out coupled in each direction with 14 bits (2^{14} levels) of precision. The background noise, such as the arbitrary reflections from the sample surface and the sensor thermal noise, was minimized using noise-canceling method and post image processing, thence the light intensity from the grating coupler regions in each direction was acquired (more details are provided in Section 4) (Figure 2C). As the result, the light intensity from all three output ports have a ratio of 22.3% : 23% : 22.1% (or 1744.5 : 1801.6 : 1727.2 from the image pixel readout) which is in excellent agreement with the FDTD simulation result with less than 0.5% deviation. Nevertheless, we envision that high-speed, low noise germanium [20], [21], [22] or graphene photodetectors [23], [24], [25], [26] can be integrated into the device and used for improving both detectability and data collection speed and accuracy.

### 2.2 Solving PDE using Kirchhoff’s photonic nodes

To showcase the functionality of a photonic platform formed by KPNs, we aim to approximate a finite difference node, which locally discretizes a Laplace equation using a finite difference method (more details in the Supplementary material, Sections 2–3). Although, this is not the only possible application space for the proposed architecture since a network comprising of KPNs could be used as recurrent neural network [27] or as a compact solution for reconfigurable routing and network broadcasting [28], or also simply as a reprogrammable filter for information pre-processing applications such as for network-edge devices.

As a proof of principle, here we select a two-dimensional heat transfer problem represented by a steady-state Laplace’s homogeneous equation (Figure 3A.i), which can be mathematically described by Eq. (1), which describes the relation between a variable *f* and its partial derivatives. Typically, PDEs are solved numerically by discretizing space (and/or time) into meshes points, in such a way that the partial derivatives can be reduced into linear combinations of the variable values at several neighboring nodes of the mesh.

In details, after applying the Finite Difference Method (FDM) to a mesh network (Figure 3A.ii), the central node *O*_{i,j} can be represented by its four adjacent nodes (Eq. (2)), where *h*_{i} is the mesh step that describes the discretization level of the problem in the analytical domain. Once the discretized mesh node is set with node-to-node correlation function *є*_{i} approximate to a constant value when the equidistant mesh step *h* is small enough, this Laplace’s equation can be locally converted in summation of incremental ratios of physical quantities and solved iteratively (further details in Supplementary material, Section 2).

However, this usually requires a large amount of compute power, memory, and scales exponentially as the problem size and required accuracy.

As suggested more than six decades ago [1], electrical resistor networks (Figure 3A.iii) can one-to-one map finite difference mesh grids. Similar to the analytical model, electrical current *I*_{i} and resistance *R*_{i} can be mapped to such a node model ruled by Kirchhoff’s law and Ohm’s law written as *I*_{i} = (*V*_{i} − *V*_{o})/*R*_{i}. Here, the current injected into the selected node always equals to the current leaving the node (Kirchhoff’s law), while the current splitting ratios to each direction will be automatically “adjusted” by the intrinsic electrical potentials of each path according to current (voltage) partitions, thus providing a solution to the problem (e.g., PDE) (more details in Supplementary material, Sections 2–3).

While the mapping of the finite difference algorithm in an electrical circuit is completely homomorphic, in photonics due to its distributed nature this is not possible. In electric circuits, due to quasi-static approximation, a variation of just one resistance, which represents the node-to-node correlation in the discretized domain, induces a redistribution of the potential in the entire network (nonlocal effects), which allows to model different gradient effects in the distribution of physical quantities. In photonics instead, the optical power when flowing in the photonic mesh is not affected by perturbation in other paths in a global sense, though this effect can be forced or induced, as discussed later. However, addressing all the splitting ratios at each node could be unpractical under certain circumstances, adjusting only the “key nodes” (e.g., nodes on the boundaries or high-loss nodes) which located adjacent to the places that the node connectivity has relatively high variations would already bring the accuracy to an acceptable range (a boundary-weighted example and the method of generating key node configuration library is given in Sections 6 and 7 of the Supplementary material in more details). In this view, a photonic processor, with induced homomorphism, could generate an initial, low-precision and approximate solution, to be used then by an integrated high-precision digital solver which ultimately produce the required high-precision solutions. This type of photonic engine can, indeed, be configured as an accelerator as it would reduce the number of iterations required by the iterative solver when solving partial differential equations associated to complex numerical simulation, e.g., Newton’s method. Practically, for this kind of solvers, it is essential to obtain an initial approximation to the solution which can be used as an initialization stage for obtaining a simulation that converges much faster.

Nevertheless, besides the disadvantage of populating the space of configurations to mimic different lumped circuit behaviors, using integrated photonics would provide three immediate advantages over resistive networks. First, the absence of charge/discharge of the wires enables distributed communication and concurrent low power dissipation. Second, once the network is set, the picosecond-short delay is dominated by the time-of-flight of the photon within the PIC. Third, the amplitude (power) of the traveling wave in the waveguides can be easily tuned using attojoule-efficient modulators through wisely engineered light–matter interactions [29], [30], [31], [32] mimicking different configurations (programmability). For the sake of simplicity, we demonstrate that it is possible to obtain accurate solutions for a uniform domain, but this approach could be extended exploiting the network reconfigurability. In the supporting information, exploiting the symmetry and reciprocity of the KPN, we show that the same chip can be used for obtaining the solution of the same Laplace PDE with different Dirichlet’s boundary conditions (Supplementary material, Section 5, Figure S9).

After providing a practical exhibition and guidelines for obtaining a Kirchhoff’s equivalent law in photonics and consequently obtaining an FDM-like node, we cascaded multiple nodes building a 5 × 5 optical FDM mesh grid to solve a discretized heat transfer problem. The assembled system maps a symmetric type of heat transfer problem with a heat source injected from the center-left of the mesh grid and surrounded by constant temperature boundary.

The input signal, which in this case represents the Dirichlet’s boundary conditions may, in general, be any arbitrary laser beam distribution coupled into any node of the circuit (i.e., here grating couplers used, Figure 3B). For mimicking a 3 × 3 FDM mesh a 5 × 5 photonic mesh is fabricated and tested, in which the additional nodes on the sides of the domain are used for reading/applying the optical power at the boundary. In electrical circuits, a Dirichlet’s boundary condition is provided by a constant potential which could be either “active” (heat source), using electrical sources (voltage or current generators), or “passive” (heat sink) as electrical ground. The electrical node at the boundary has two functions: (1) sets the value of the function at the boundary by applying a constant voltage (2) due to the lumped nature of the circuit, forces a gradient trend (differential voltage) throughout the lumped circuit. In a photonic network, if the former is possible by applying a fixed optical power, the latter is not straightforwardly achievable due to the absence of nonlocal effects. In this view, this represents another breaking point in the likeness between electrical and optical meshes. Active boundary conditions, mimicked by local sources, sets the value of the optical power at boundary, while passive boundary conditions (no laser source applied) do not force the optical power at the boundary node to be zero. This is achieved in the successive normalization step. Additionally, the gradient trend, which in the electrical mesh is stringently related to the boundary conditions, in SPACE is fully decoupled. Our approximation consists in directly embedding into the network the boundary-induced thermal distribution by modifying the node-to-node optical losses (i.e., adding 1 dB of loss by using waveguide bending, like in this case, or actively as shown in Section 6 of the Supplementary materials), accounting for modified node-to-node correlation function *є*_{i}. This is a further step towards our quest to an induced homomorphism.

To characterize the performance of the system and obtaining discretized measurements for each node, first, we introduce for each direction of the nodes a set of 50/50 Y-branch splitter followed by a grating coupler in order to estimate the optical power at each node. The power drop at each node represents the temperature distribution at each point of the discretized domain, and it is measured, as previously observed, in time-parallel through a properly calibrated camera (Figure 3C). However, for high-speed reconfiguration operation, integrated photodetectors can be used with a latency on the order of 10’s ps.

In order to obtain readable data from the furthest node from the input, 39 mW of laser input (as the maximum power output from our laser source) is applied to the 5 × 5 SPACE mesh grid. The optical nonlinearity effect has also been considered and the actual optical power coupled into the first node is well below 5 mW to prevent this [33].

We verify the accuracy of the approximate solution of the 5 × 5 SPACE prototype by comparing the obtained experimental measurement at each KPN to the discretized and normalized solution of a heat transfer problem obtained through ac commercially available numerical solver (COMSOL Multiphysics) with the same mesh resolution and density. Considering the different input units (i.e., temperature and optical power) and values (i.e., 100 °C in the thermal model and 39 mW optical power in the measurement), all the output at each node are normalized with respect to the input, thus allowing a dimensionless comparison.

A discretized solution is obtained through commercially available software (COMSOL) which served as a problem accuracy baseline. It is worth mentioning that any discretized solution is proportional to the mesh resolution (e.g., a 5 × 5 mesh COMSOL simulation has 99.95% accuracy comparing to a 300 × 300 mesh averaged down to a 5 × 5 with same initial setup, Figure S3). In details, in order to decouple the discretization error from the error produced by the approximation of the model, a 5 × 5 discretized mesh is used as comparison with the experimentally measured solutions at each node of the 5 × 5 SPACE (Figure 4A). SPACE provides an overall average accuracy of 98% and a discrete solution with an accuracy higher than 95% for each node (Figure 4B).

In our pursuit of a forced homomorphism, similarly to the tolerance of the resistors in a resistive network that mimics a FDM [1], the deviations caused by fabrication variances, such as the grating coupler efficiency, *y*-branch splitting ratio, and variability of coupling coefficients perceptively lower the accuracy of the solution.

Due to the distributed nature of the photonic engine, several iterations, consisting of multiple reflections at each discontinuity, are needed in SPACE to obtain a time-stable solution. This time is proportional to the size of the network and the density, i.e., number of KPN. Based on PIC simulations performed using Lumerical Interconnect, a full iteration cycle (dominated by the time of flight of the photon) takes only 30 ps considering 100 μm node-to-node (*n*2*n*) spacing for a 5 × 5 SPACE in this sparse design. The iteration time drops to 16 ps with 25 μm *n*2*n* spacing which is the highest density achievable in our current design (Figure 4C). The highest density can be achieved by using on-chip integrated photodetectors as a detection mechanism, instead of a IR-camera which requires outcoupled radiation from the chip by means of taps and grating couplers at each node which compel sufficient space allocated on chip. Furthermore, considering the photonic node size, the SPACE engine can be packed with a minimum density of 25 μm/component, although we separate the nodes with 200 μm spacing for reducing the output crosstalk while measuring. The footprint of the network can be further shrunk using inverse design approaches obtaining a density of <5 µm/node, enabling higher density meshes. Additionally, SPACE can be fabricated with adaptive mesh, with an increased density within certain sensitive or turbulent regions of the simulation, thus increasing the overall accuracy. It is also worth to mention that, as an approximate computing engine, when the target accuracy is relaxed to 90% of its maximum, the iteration time drops to 1.8 ps which is equivalent to 556 GHz. In terms of the scalability, the latency does increase as the network size scales from 5 × 5 to 10 × 10, in which the light propagation time in the waveguide will contribute even less to the total runtime, thus proving that SPACE could be potentially further scaled-up (further details in Supplementary material, Section 5).

Our KPN design and SPACE circuit provide a powerful tool for homogenously distributing optical power in a defined network similarly to a lumped circuit subjected to Kirchhoff’s Law. When the correlation function between each node of the network is wisely selected or actively tuned (e.g., with electro-optic modulators or switches), the network, mimicking FDM, can solve a general second-order Laplace’s PDE in an analog manner. For instance, if the splitting ratio of nodes adjacent to the boundary conditions could be tuned to 10% : 10% : 80% (i.e., more light routes to the boundaries), a 5 × 5 modulated SPACE is able to improve the accuracy up to 99.2% (more details discussed in Supplementary material, Section 6). Similar configurations may potentially be explored to solve time dependent or nonlinear PDE with arbitrary boundary conditions by introducing time discretization or nonlinear elements, respectively. PDEs applied to nonhomogeneous domains can be mapped on SPACE by changing the splitting ratio according to characteristic distribution (e.g., different thermal conductivity mapped as different attenuation for each individual node). Other cases of PDE applied to nonsymmetrical or inhomogeneous domains are reported as PIC numerical simulation in the Supplementary material.

Towards dynamic problem reconfigurability, if electro-absorption modulators are introduced between neighboring nodes in SPACE, a vast number of different PDEs can be solved. For illustrative purposes, we numerically show, using a photonic interconnection emulator, that it is possible to map the temperature distribution, solution of Laplace equation, onto SPACE by adjusting the extinction ratio of these modulators between neighboring nodes according to the problem to be solved. This allows solving a multitude of problem cases for number of exemplary underlying heat-conducting materials (details in the Supplementary material, Section S7). Similar to commercial numerical solvers, we present how KPN and SPACE with added reconfigurability produces look-up table solutions for the specific configuration.

Beyond the exemplary Laplace equation investigated thus far, other PDEs can, in principle, be solved such as Poisson’s equation, for example, if additional light sources are added to the nodes, mimicking the different node potential [34]. Nonetheless, other PDEs like diffusion equations and wave equations would require optical capacitive and inductive elements needed to express the time-dependent variances enabling the one-shot solution. Different from emulating an optical resistor, which can be easily realized by optical lossy materials or electro-optical modulators, optical capacitors and inductors require specific designs to mimic the behavior of their electrical counterpart. For example, a Fabry–Perot interferometer with chirped Bragg gratings have been demonstrated as an optical capacitive component which can act like a broadband low pass or high pass filter [35]. On the other hand, an optical inductive component can be implemented as a self-electro-optical device with both integrated modulator and detectors that use the photocurrent to back feed the modulator and change the light intensity injected into the detection region as a negative feedback loop [36], [37]. Reconfigurability of the optical mesh nodes, and hence boundary condition programmability, can also be introduced via electro-optic tunable materials such the unity-strong index variation demonstrated in transparent conductive oxides (e.g., ITO [9], [30]) with temporal speeds above 1 GHz [38]. Indeed, control of sub-wavelength small programmable elements allows to realize optical circuitry enabling designing PDE solvers that are actually governed by the quasi-static approximation of Maxwell’s equation, thus by Kirchhoff’s law [39]. In such optical circuits, the effective wavelength is “stretched” due to the effect of epsilon near (or at) zero media leading to nonlocal interactions within the optical circuit board. The propagating (conduction) electric displacement inside such optical nano-circuits can then be probed and provide, in principle, error-free, solutions to the PDE problem. Realistic implementations of such nano-optical analog processors, however, need to address the nonvanishing imaginary part of the optical index, if higher than about 95% solution accuracy is desired [39].

## 3 Conclusions

In summary, we propose the designs of a photonic node which is able to replicate the equivalent of Kirchhoff’s law for optical power. Using equal splitting functionality, we replicate a mesh structure which approximate a homogeneous lumped circuit model. We use the photonic circuit to map a finite difference approach to solve partial differential equation effortlessly and noniteratively, termed Silicon Photonic Approximate Computing Engine (SPACE). Our numerical and experimental analysis indicates that the steady-state response of a characteristic SPACE engine may be achievable in 16 ps, obtaining inherently discretized solutions for each point of the mesh of the domain with a bandwidth up to 63 GHz. This approach could easily adapt another active component, such as electro-optic modulators, photodetectors, and tunable photonic cavities like photonic crystals, to solve more complicated problems with denser, heterogeneous meshes and arbitrary boundary conditions.

Furthermore, the proposed approach features reconfigurability of the input positions and boundary conditions with low-loss network interconnectivity of such distributed networks via PICs and ensures foundry-near cost scaling. This approach allows solving more intricate problems such as those with heterogeneous grid and Neumann boundary conditions when chip is augmented by active components, such as electro-optic modulators, photodetectors, and tunable photonic cavities like photonic crystals. Our findings provided a novel pathway to ultra-fast, integrable, and reconfigurable photonic analog computing engines based on an integrated photonics Kirchhoff’s node, used for solving PDEs, but can also be adopted as core structure in recurrent neural networks or as compact solution for network broadcasting.

## 4 Methods

### 4.1 Fabrication process

All of the single optical power splitters and 5 × 5 SPACE is fabricated on the same 220 nm Silicon on Insulator (SOI) chip to minimize the variance during the fabrication process. Raith Voyager 50 kV E-beam lithography system is used with fix beam moving stage (FBMS) feature to allow zero waveguide stitching errors across multiple write fields. Hydrogen silsesquioxane (HSQ) with 2% concentration is used to provide around 42 nm of mask thickness (4000 rpm for 60 s) with high resolution in writing. After the spin coating, the chip is put on a hotplate for 240 s pre-bake at 80-degree centigrade. After the patterning, the chip is dipped into MF-319 for 70 s to develop the unexposed HSQ area including 5 s of gentle stirring to shake off the air bubbles of the chemical reaction. Then 30 s of D.I. water rinse will be immediately applied to stop the development and clean up the residue. To etch down the silicon layer and reveal the features, a 28 s of SF_{6} and C_{4}F_{8} (both at 10 sccm) at 500 W ICP power and 20 W bias etching with Plasma-Therm Apex SLR Inductively Coupled Plasma Etcher is able to fully etch all the silicon down and provide over 9:1 selectivity for our smallest features.

### 4.2 Measurement and data processing

To measure the output light intensity, an optical probe station setup is used with a tunable laser at 1550 nm wavelength connecting to a lens fiber to maximize the light coupled onto the chip. Considering the polarization of the grating coupler and its coupling efficiency, the actual laser power coupled into the mesh is less than 5 mW, which is still far below the nonlinearity energy density limitation of the Silicon Photonic waveguide (500 nm × 220 nm). Xenics IR camera integrated with the microscope captures the scattering light at each output grating. In addition, a black light shield is applied to cover the entire camera, probe station and microscope to prevent the ambient light. And the thermal noise of the camera is eliminated by capturing the image with no laser input. The last type of noise taken into account in the measurement is the surface reflection including the lens flare, and this is by substituting the averaged background readout that adjacent to the grating coupler. After the noise cancelation, the images are imported into MATLAB to integrate the intensity values (0–4095 for our 12-bit depth sensor) of all the pixels of the output region. It is also worth to mention that nodes at different positions have over three orders of magnitude difference which is far beyond the dynamic range of the camera. Therefore, lower input laser power with shorter camera integration time is used for nodes closer to the input node and post-processed into the same scale.

**Funding source: **National Science Foundation

**Award Identifier / Grant number: **1748294

**Author contributions:**V.J.S. and T.E. conceived the original ideas, acquired the funds, and supervised the project. S.S. and M.M. developed the relevant theories and analyses for the project. S.S. and M.M. designed the nano-optic circuit. S.S. fabricated and tested the optical engine and S.S., X.M. and Z.M. conducted the related numerical simulations. C.S. simulated the reprogrammable photonic engine, while E.K. and J.A. evaluated and benchmark the performance of the photonic chip. M.M., S.S., T.E. and V.J.S. discussed the results and contributed to the understanding, analysis, and interpretation of the results. MM wrote the first draft of the manuscript, and V.J.S., and S.S. contributed to writing subsequent drafts of the manuscript. All authors provided critical feedback.**Research funding**: This work was supported by National Science Foundation under award number (1748294).**Conflict of interest statement:**The authors have filed a provisional patent application on this idea, U.S. Patent 10,318,680. No other conflict of interest is present.**Data and materials availability:**All data needed to evaluate the conclusions in the paper are present in the paper or the supplementary materials.

## References

[1] G. Liebmann, “Solution of partial differential equations with a resistance network analogue,” *Br. J. Appl. Phys.*, vol. 1, pp. 92–103, 1950, https://doi.org/10.1088/0508-3443/1/4/303.Search in Google Scholar

[2] D. C. Barker, “Electrical analogue for a partial differential equation,” *Math. Comput. Simulat.*, vol. 16, pp. 38–45, 1974, https://doi.org/10.1016/s0378-4754(74)80013-3.Search in Google Scholar

[3] B. S. Vimoke, *Simulating Water Flow in Soil with an Electrical Resistance Network*, USA, Agricultural Research Service, U.S. Department of Agriculture, 1964.Search in Google Scholar

[4] G. M. Amdahl, in *Proceedings of the April 18–20, 1967, Spring Joint Computer Conference. AFIPS ’67 (Spring)*, Association for Computing Machinery, Atlantic City, New Jersey, 1967, pp. 483–485. https://doi.org/10.1145/1465482.1465560.Search in Google Scholar

[5] A. Adovic, M. Williams, D. Rousseau, et al., “Towards systems-on-a-chip,” *Nat. Photonics*, vol. 12, 2018, Art no. 311.10.1038/s41566-018-0191-7Search in Google Scholar

[6] X. Qiang, X. Zhou, J. Wang, et al.., “Large-scale silicon quantum photonics implementing arbitrary two-qubit processing,” *Nat. Photonics*, vol. 12, pp. 534–539, 2018, https://doi.org/10.1038/s41566-018-0236-y.Search in Google Scholar

[7] D. A. B. Miller, “The role of optics in computing,” *Nat. Photonics*, vol. 4, p. 406, 2010, https://doi.org/10.1038/nphoton.2010.163.Search in Google Scholar

[8] M. Miscuglio, A. Mehrabian, Z. Hu, et al.., “All-optical nonlinear activation function for photonic neural networks [Invited],” *Opt. Mater. Express*, vol. 8, pp. 3851–3863, 2018, https://doi.org/10.1364/ome.8.003851.Search in Google Scholar

[9] R. Amin, J. K. George, S. Sun, et al.., “ITO-based electro-absorption modulator for photonic neural activation function,” *Apl. Mater.*, vol. 7, 2019, Art no. 081112. https://doi.org/10.1063/1.5109039.Search in Google Scholar

[10] Y. Shen, N. C. Harris, S. Skirlo, et al.., “Deep learning with coherent nanophotonic circuits,” *Nat. Photonics*, vol. 11, pp. 441–446, 2017, https://doi.org/10.1038/nphoton.2017.93.Search in Google Scholar

[11] A. N. Tait, T. F. de Lima, E. Zhou, et al.., “Neuromorphic photonic networks using silicon photonic weight banks,” *Sci. Rep.*, vol. 7, pp. 1–10, 2017, https://doi.org/10.1038/s41598-017-07754-z.Search in Google Scholar

[12] N. Mohammadi Estakhri, B. Edwards, and N. Engheta, “Inverse-designed metastructures that solve equations,” *Science*, vol. 363, pp. 1333–1338, 2019, https://doi.org/10.1126/science.aaw2498.Search in Google Scholar

[13] J. Wu, P. Cao, X. Hu, et al.., “Compact tunable silicon photonic differential-equation solver for general linear time-invariant systems,” *Opt. Express*, vol. 22, pp. 26254–26264, 2014, https://doi.org/10.1364/oe.22.026254.Search in Google Scholar

[14] S. Tan, L. Xiang, J. Zou, et al.., “High-order all-optical differential equation solver based on microring resonators,” *Opt. Lett.*, vol. 38, pp. 3735–3738, 2013, https://doi.org/10.1364/ol.38.003735.Search in Google Scholar

[15] J. Anderson, E. Kayraklioglu, S. Sun, et al., “ROC: A reconfigurable optical computer for simulating physical processes,” *ACM Trans. Parallel Comput.*, vol. 7, no. 1, 2020. https://doi.org/10.1145/3380944.Search in Google Scholar

[16] V. J. Sorger, S. Sun, T. El-Ghazawi, A.-H. A. Badawy, and V. K. Narayana, *Reconfigurable Optical Computer*, 2017. Available at: https://patents.google.com/patent/US20170161417A1/en.Search in Google Scholar

[17] J. Lu and J. Vučković, “Nanophotonic computational design,” *Opt. Express.*, vol. 21, 2013, Art no. 13351. https://doi.org/10.1364/oe.21.013351.Search in Google Scholar

[18] A. Y. Piggott, J. Lu, K. G. Lagoudakis, J. Petykiewicz, T. M. Babinec, and J. Vučković, “Inverse design and demonstration of a compact and broadband on-chip wavelength demultiplexer,” *Nat. Photonics*, vol. 9, pp. 374–377, 2015, https://doi.org/10.1038/nphoton.2015.69.Search in Google Scholar

[19] S. Molesky, Z. Lin, A. Y. Piggott, W. Jin, J. Vucković, and A. W. Rodriguez, “Inverse design in nanophotonics,” *Nat. Photonics*, vol. 12, pp. 659–670, 2018, https://doi.org/10.1038/s41566-018-0246-9.Search in Google Scholar

[20] H. Chen, P. Verheyen, P. D. Heyn, G. Lepage, J. D. Coster, P. Absil, G. Roelkens, and J. V. Campenhout, in *2016 Optical Fiber Communications Conference and Exhibition (OFC)*, 2016, pp. 1–3.Search in Google Scholar

[21] L. Chen, P. Dong, and M. Lipson, “High performance germanium photodetectors integrated on submicron silicon waveguides by low temperature wafer bonding,” *Opt. Express*, vol. 16, 2008, Art no. 11513. https://doi.org/10.1364/oe.16.011513.Search in Google Scholar

[22] *High-Performance Ge-on-Si Photodetectors | Nature Photonics*. Available at: https://www.nature.com/articles/nphoton.2010.157.Search in Google Scholar

[23] X. Gan, R.-J. Shiue, Y. Gao, et al., “Chip-integrated ultrafast graphene photodetector with high responsivity,” *Nat. Photonics*, vol. 7, 2013, Art no. 883. https://doi.org/10.1038/nphoton.2013.253.Search in Google Scholar

[24] I. Goykhman, U. Sassi, B. Desiatov, et al., “On-chip integrated, silicon–graphene plasmonic Schottky photodetector with high responsivity and avalanche photogain,” *Nano Lett.*, vol. 16, pp. 3005–3013, 2016, https://doi.org/10.1021/acs.nanolett.5b05216.Search in Google Scholar

[25] A. Pospischil, M. Humer, M. M. Furchi, et al., “CMOS-compatible graphene photodetector covering all optical communication bands,” *Nat. Photonics*, vol. 7, 2013, Art no. 892. https://doi.org/10.1038/nphoton.2013.240.Search in Google Scholar

[26] Z. Ma, K. Kikunaga, H. Wang, et al., “Compact graphene plasmonic slot photodetector on silicon-on-insulator with high responsivity,” *ACS Photonics*, 2020, https://doi.org/10.1021/acsphotonics.9b01452.Search in Google Scholar

[27] T. W. Hughes, I. A. D. Williamson, M. Minkov, and S. Fan, “Wave physics as an analog recurrent neural network,” *Sci. Adv.*, vol. 5, 2019, Art no. eaay6946. https://doi.org/10.1126/sciadv.aay6946.Search in Google Scholar

[28] V. K. Narayana, S. Sun, A.-H. A. Badawy, V. J. Sorger, and T. El-Ghazawi, “MorphoNoC: exploring the design space of a configurable hybrid NoC using nanophotonics,” *Microprocess. Microsyst.*, vol. 50, pp. 113–126, 2017, https://doi.org/10.1016/j.micpro.2017.03.006.Search in Google Scholar

[29] D. A. B. Miller, “Attojoule optoelectronics for low-energy information processing and communications,” *J. Lightwave Technol.*, vol. 35, pp. 346–396, 2017, https://doi.org/10.1109/jlt.2017.2647779.Search in Google Scholar

[30] R. Amin, R. Maiti, C. Carfano, et al., “52 V mm ITO-based Mach–Zehnder modulator in silicon photonics,” *APL Photonics*, vol. 3, 2018, Art no. 126104. https://doi.org/10.1063/1.5052635.Search in Google Scholar

[31] R. Amin, Z. Ma, R. Maiti, et al., “Attojoule-efficient graphene optical modulators,” *Appl. Opt.*, vol. 57, 2018, Art no. D130. https://doi.org/10.1364/ao.57.00d130.Search in Google Scholar

[32] V. J. Sorger, R. Amin, J. B. Khurgin, Z. Ma, H. Dalir, and S. Khan, “Scaling vectors of attoJoule per bit modulators,” *J. Opt.*, vol. 20, 2018, Art no. 014012. https://doi.org/10.1088/2040-8986/aa9e11.Search in Google Scholar

[33] H. K. Tsang and Y. Liu, “Nonlinear optical properties of silicon waveguides,” *Semicond. Sci. Technol.*, vol. 23, 2008, Art no. 064007. https://doi.org/10.1088/0268-1242/23/6/064007.Search in Google Scholar

[34] G. B. Folland, *Introduction to Partial Differential Equations*, 2nd ed. Princeton, N.J., Princeton University Press, 1995.Search in Google Scholar

[35] S. Loranger, M. Gagné, and R. Kashyap, “Capacitors go optical: wavelength independent broadband mode cavity,” *Opt. Express*, vol. 22, pp. 14253–14262, 2014, https://doi.org/10.1364/oe.22.014253.Search in Google Scholar

[36] A. Majumdar and A. Rundquist, “Cavity-enabled self-electro-optic bistability in silicon photonics,” *Opt. Lett.*, vol. 39, pp. 3864–3867, 2014, https://doi.org/10.1364/ol.39.003864.Search in Google Scholar

[37] S. Sun, R. Zhang, J. Peng, et al., “MO detector (MOD): A dual-function optical modulator-detector for on-chip communication,” *Opt. Express*, vol. 26, 2018, Art no. 8252. https://doi.org/10.1364/oe.26.008252.Search in Google Scholar

[38] R. Amin, R. Maiti, Y. Gui, et al., “Sub-wavelength GHz-fast broadband ITO Mach–Zehnder modulator on silicon photonics,” *OPTICA*, vol. 7, pp. 1812–1819, 2020, https://doi.org/10.1364/optica.389437.Search in Google Scholar

[39] M. Miscuglio, Y. Gui, X. Ma, et al., *Analog Computing with Metatronic Circuits*, arXiv preprint: 2007.05380.Search in Google Scholar

## Supplementary Material

The online version of this article offers supplementary material (https://doi.org/10.1515/nanoph-2020-0655).

**Received:**2020-12-15

**Accepted:**2021-03-02

**Published Online:**2021-03-22

© 2021 Shuai Sun et al., published by De Gruyter, Berlin/Boston

This work is licensed under the Creative Commons Attribution 4.0 International License.