High-throughput optical neural networks based on temporal computing
An emerging generative artificial intelligence (AI) based on neural networks
starts to grow in popularity with a revolutionizing capability of creating new
and original content. As giant generative models with millions to billions of
parameters are developed, trained and maintained, a massive and
energy-efficient computing power is highly required. However, conventional
digital computers are struggling to keep up with the pace of the generative
model improvements. In this paper, we propose and demonstrate high-throughput
optical neural networks based on temporal computing. The core weighted
summation operation is realized with the use of high-speed electro-optic
modulation and low-speed balanced photodetection. The input data and weight are
encoded in a time sequence separately and loaded on an optical signal via two
electro-optic modulators sequentially. By precisely controlling the
synchronization time of the data and weight loading, the matrix multiplication
is performed. Followed by a balanced photodetector, the summation is conducted,
thanks to the electron accumulation of the inherent electronic integrator
circuit of the low-speed photodetector. Thus, the linear weighted summation
operation is implemented based on temporal computing in the optical domain.
With the proposed optical linear weighted summation, a fully-connected neural
network and convolutional neural network are realized. Thanks to the high-speed
feature of temporal computing, a high data throughput of the optical neural
network is experimentally demonstrated, and the weighting coefficients can be
specified on demand, which enables a strong programmability of the optical
neural network. By leveraging wavelength multiplexing technology, a scalable
optical neural network could be created with a massive computing power and
strong reconfigurability, which holds great potential for future giant AI
applications.
Shuang Zheng, Jiawei Zhang, Weifeng Zhang
## Scalable optical neural networks based on temporal computing

## Abstract

The optical neural network (ONN) has been considered as a promising candidate for next-generation neurocomputing due to its high parallelism, low latency, and low energy consumption, with significant potential to release unprecedented computational capability. Large-scale ONNs could process more neural information and improve the prediction performance. However, previous ONN architectures based on matrix multiplication are difficult to scale up due to manufacturing limitations, resulting in limited scalability and small input data volumes. To address this challenge, we propose a compact and scalable photonic computing architecture based on temporal photoelectric multiplication and accumulate (MAC) operations, allowing direct processing of large-scale matrix computations in the time domain. By employing a temporal computing unit composed of cascaded modulators and time-integrator, we conduct a series of proof-of-principle experiments including image edge detection, optical neural networks-based recognition tasks, and sliding-window method-based multi-target detection. Thanks to its intrinsic scalability, the demonstrated photonic computing architecture could be easily integrated on a single chip toward large-scale photonic neural networks with ultrahigh computation throughputs.
## Introduction

With the advent of the big-data era, the proliferation of artificial intelligence (AI) applications and powerful machine learning techniques create exponentially increased data streams that can hardly be processed efficiently by conventional computing systems and architectures[1, 2]. In the meanwhile, the unsustainability of Moore's law and the failure of Dennard Scaling rules greatly hinder the development of current computing hardware based on electronic transistor and von Neumann architecture[3, 4]. To resolve this disagreement, optical computing has been envisioned as a promising solution, which holds great potential to release unprecedented computational capability by leveraging the advantages inherent to optical signals, such as low power consumption, low latency, ultra-wide bandwidth and high parallelism[5, 6, 7]. In particular, in past years, optical neural network (ONN) layouts have already penetrated successfully the deep learning era. With the rapidly increasing complexity in dataset size of various application scenarios, large-scale ONN is highly desired to meet technological demands and performance criteria. To achieve better accuracy and efficiency in complex neural tasks, it is critical to scale up the dimensions of ONNs, including depth, width, and resolution, which correspond to the number of layers, the number of neurons in a layer, and the size of the input image, respectively[8, 9]. Due to the advantages of large-scale matrix operation and high computational capability, spatial optical computation has shown great potential for establishing large and compact computing units predominantly required in optical-artificial intelligence computers[10, 11, 12]. Fourier transforms and convolutions as well as applications in classification tasks have been proposed utilizing the optical diffraction masks in free space. However, traditional free-space approaches require heavy auxiliary or modulation equipment[13, 14, 15, 16, 17], such as lenses, spatial light modulators (SLMs) and digital micro-mirror devices (DMDs), which may make the system bulky and costly, hindering their feasibilities towards the large-scale implementation and commercialization of photonic neural networks. Therefore, there remain great opportunities for substantial improvements in photonic hardware accelerators.

With the rapid development of advanced photonic manufacturing technique, silicon photonic integrated circuits (PICs) have been considered as a promising candidate for establishing large-scale photonic neural networks, which have the advantages of compact size, high-density integration, and power efficiency[18, 19, 20, 21, 22]. Typical silicon PIC-based ONN architectures use cascades of multiple Mach-Zehnder interferometers (MZIs) and microring resonator array-based wavelength division multiplexing (WDM) technology. The number of units for these ONN architectures scales quadratically with the input data dimension for matrix multiplication, leading to quadratic scaling in the footprint and energy consumption. Though PICs have shown great potential in integrated optics, the integration scale due to manufacturing limitations, and the complex control circuits for reducing fluctuation of the resonance wavelength restrict the development of a large-scale programmable photonic neural network. Consequently, the size of the input vectors and the
matrices that can be loaded on these hardwares is extremely limited by the integration scale. More recently, chip-scale microcomb devices have been employed as multi-wavelength light sources [37; 38; 39; 40], which can bring much convenience and greatly increase the input data size. However, the input matrix size may be still limited by the number of multiplexed wavelength channels. Additionally, there are also other on-chip architectures that can support matrix loading and computing [41; 42; 43; 44; 45]. Despite impressive performance, their practical compute rate and processed matrix size are still far beyond their promise for high-speed photonic engines. In this scenario, an integrable and scalable optical computing architecture featuring large-scale matrix operation is highly desired, which would increase the network size and resolution to enhance predictive accuracy.

In this paper, we propose an integrable and scalable photonic computing architecture capable of performing parallel MAC operations. Different from prior works, the MAC operations are implemented in the time domain by using a fundamental execution unit. The intrinsic temporal accumulation mechanism of the proposed approach determines that the temporal, wavelength and spatial dimensions can be fully utilized to perform large-scale matrix computation and enhance the computation capability [46; 47]. We conduct a series of proof of concept experiments to take full advantage of the proposed computing architecture, including typical image edge detection, fully-connected neural network and convolutional neural network-based classification tasks, and a primary validation for sliding-window method-based multi-target detection task. With further on-chip integration and multidimensional multiplexing, the demonstrated photonic computing architecture is promising toward large-scale photonic neural networks with ultrahigh computation throughputs.
## Results

**Scalable photonic computing architecture.** The key components in the proposed scalable photonic computing archihtecture are the temporal photonic MAC operation units, which are used as the fundamental

Figure 1: **Scalable photonic computing architecture.** (a) Schematic of the proposed scalable photonic computing architecture, which is implemented by employing temporal, wavelength and spatial dimensions. When the light beams transmit through cascaded MZMs, the multiplication process is accomplished. The temporal waveforms are detected by BPDs and time-integrators, effectively achieving the temporal accumulation operation. (b) For the input data, the matrix rows and columns are mapped to a time-wavelength basis. The input row vectors are encoded as temporal serial symbols and loaded onto different light wavelengths by high-speed modulators. The weight matrix rows and columns are mapped to a space-time basis. The column vectors are also encoded as temporal serial symbols by another MZM array. (c) Schematic of single-channel temporal photonic MAC operation unit. In the proposed solution, the input matrix could be encoded as non-negative values, and the weight matrix could be encoded as both positive and negative values. LD, laser diode; MUX, wavelength division multiplexer; BS, beam splitter; DEMUX, wavelength division demultiplexer; VOA, variable optical attenuator, BPD, balanced photodetector.
computing units to multiply and accumulate temporal signals. Figure 1(a) shows the proposed scalable photonic computing architecture, which could perform large-scale matrix operations by simultaneously exploiting temporal, wavelength and spatial dimensions. The architectures mainly consists of three major parts. The first part is composed of an MZM array, which is used for high-speed image data loading. The data matrix rows and columns are mapped to a time-wavelength basis. As shown in Fig. 1(b), the input row vectors of size N are encoded as temporal serial symbols and loaded onto M different light wavelengths by high-speed modulators, thus an M\(\times\)N data matrix with M rows and N columns could be loaded. All the wavelengths are then combined through a wavelength division multiplexer, and fanned out to L different spatial channels. Thus, M\(\times\)L multiplexed channels can be fully used to enhance both the weight matrix dimension and the computational capability. In the second part, another MZM array is used to load the weight matrix and also functions as multipliers. The weight matrix rows and columns are mapped to a space-time basis. Due to the broadband modulation characteristic of the MZMs, the values in data matrix (M\(\times\)N) could be correspondingly multiplied with the values in the column vector (N\(\times\)1) loaded in each spatial channel. The third part of the computing architecture is used for the temporal accumulation process. The multiplexed wavelengths in each spatial channel are separated into different optical paths via a wavelength division demultiplexer. Then, the temporal optical signals are converted into electrical signals by a high-speed BPD array. At last, time-integrators are used to complete the accumulation of temporal serial signals. Thanks to its temporal accumulation mechanism of the proposed method, the length of the loaded temporal vector is theoretically unlimited. At the same time, by multiplexing wavelength and spatial dimensions, large-scale matrix-vector and matrix-matrix multiplication could be realized on the proposed temporal photonic computing architecture.

Figure 1(c) shows the details of the proposed temporal photonic MAC operation unit, where the values in the weight matrix could be positive or negative. The input data is first flattened into a vector M and encoded as a time-domain serial electrical signal, and then loaded onto the optical carrier by the first-stage MZM with a half-wave voltage of \(V_{\pi}\). Thus, the input vector values are represented by the modulated light intensity \(I_{\text{M}}\), which can be expressed as \(\dfrac{M(t)}{V_{\pi}}\) after linearization. After the first-stage MZM, the optical signal \(I_{\text{M}}\) is divided into two parts through a 50:50 beam splitter. One part is sent to the second-stage MZM, and then detected by the positive port of the BPD. The other part is sent to the negative port of the BPD after adjusting the intensity through a variable optical attenuator (VOA). The second-stage MZM is used to encode the weight matrix. Since the weight matrices usually have both positive and negative values, the second-stage MZMs are biased at the quadrature transmission point and work in the linear region. The weight matrix W is also flattened into a vector and loaded by the second-stage MZM. After propagating through two-stage MZMs, the light intensity in the upper optical path can be expressed as:
\[I_{\rm up}(t)=\frac{1}{4}\,I_{\rm M}(t)\times\left[1+\frac{W(t)}{V_{\pi}}\right] \tag{1}\]

\(I_{\rm M}\) is the light intensity modulated by the first MZM, which represents the temporal signal sequence of the input data stream. \(W(t)\) is the temporal signal sequence representing weigth data. The signal \(I_{\rm up}\) modulated by two-stage MZMs is then sent to the positive port of the BPD. After propagating through the first-stage MZM and the VOA, the optical signal in another optical path could be expressed as:

\[I_{\rm down}(t)=\frac{\alpha}{2}\,I_{\rm M}(t) \tag{2}\]

The temporal signal \(I_{\rm down}\) is detected by the negative port of the BPD. Thus, the output of the BPD can be given by:

\[I_{\rm BPD}(t)=\frac{1}{2}\,\Re\times I_{\rm M}(t)\times\left[\frac{W(t)}{V_{ \pi}}+\frac{1}{2}-\alpha\right] \tag{3}\]

where \(\Re\) is the responsivity of the BPD. To get the multiplication result from the BPD, the attenuation coefficient \(a\) of the VOA should be equal to 1/2. To further achieve the temporal accumulation, the temporal serial signals are sent to a time-integrator. When the duration time of the optical signal that needs to be accumulated is less than the maximum integration time, the electrical output could be derived:

\[V_{\rm out}\,\propto\int_{\delta t}\left[\frac{1}{2}\,\Re\times M(t)\times W (t)\right]dt \tag{4}\]

where the output voltage \(V_{\rm out}\) is proportional to the result of the vector dot product.

In the experiment, the temporal accumulation process is achieved by using a low-speed BPD as an integrator (see Supplementary Note 1). A series of square light pulses with varied durations are detected by a BPD for the test of temporal integration. Since the duration of the square light pulses is less than the respond time of the BPD, the amplitude of the output electrical signal increases linearly with the pulse width. The details of the accumulation process are further illustrated in the experiment (see Supplementary Note 1). Hence, the temporal MAC operation could be accomplished by the proposed unit.
**Convolutional image processing with photonic MAC operation unit.** Figure 2(a) shows the experimental setup of the proposed temporal photonic MAC operation unit for convolutional image processing. In the experiment, we perform two-dimensional (2D) convolution on an 8-bit 92\(\times\)92 flower image to detect its edge feature. One channel of the RGB three-channel flower picture shown in Fig. 2(d) is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig. 2(a). The 92\(\times\) flower image is the same as the one in Fig. 2(a). The 92\(\times\)92 flower image is the same as the one in Fig.
is used as the input grayscale image. A 3\(\times\)3 Laplacian kernel \(\begin{bmatrix}-1&-1&-1\\ -1&8&-1\\ -1&-1&-1\end{bmatrix}\) is shown in Fig. 2(e), which could be used to realize a two-dimensional (2D) isotropic measurement of the second-order spatial derivative of an image, and highlights regions of rapid intensity change (See Supplementary Fig. 4 for more detailed information). According to the order in which the convolution kernel moves in the convolution calculation, the grayscale image is divided into a series of 3\(\times\)3 data groups, and each data group contains 9 pixels. Then, the image is flattened into a vector by the order of the groups. A certain time interval is added between each data group to ensure that different groups will not interfere with each other when the temporal signal is accumulated. In the experiment, the time interval between each data group is 16 ns. Same as the input image, the convolution kernel also needs to be flattened into a vector, and repeated periodically with a time interval of 16 ns. After the flattening operation, the image signal and the convolution kernel signal are loaded into the photonic MAC operation unit through the cascaded MZMs using an arbitrary waveform generator (AWG). The symbol rate of the signals generated by the AWG is 2.88 gigabaud per second, with a resolution of 8 bits per symbol. When the light propagates through the cascaded MZMs, the values of the vectors that representing the image information and the convolution kernel are multiplied in the time domain. The multiplication results for each data group are then accumulated in the time domain through a BPD with a bandwidth of 150 MHz. The processed image could be recovered by sampling and rearranging the output temporal signal from the BPD. Since there is a fixed time interval between each data group in the temporal accumulation process, the effective computing speed is 2\(\times\)9/16=1.125 GMAC operations per second. Figure 2(b) shows part of the measured temporal waveforms for convolutional image processing, and Fig. 2(c) shows the details of the electrical signal waveform output from BPD. As can be seen in Fig. 2(f) and (g), the experimentally measured edge-detected image agrees very well with the theoretically calculated one using convolutional image processing algorithms. The measured experimental result confirms the feasibility of convolutional image processing using the photonic MAC operation unit, which is a prerequisite for building a photonic convolutional neural network (OCNN).

**OCNN for image recognition.** In deep learning, a CNN is a class of ANN, most commonly applied to analyze visual imagery[49]. The powerful feature extraction capability of CNN makes it highly efficient for various image recognition tasks. In general, a typical CNN is comprised of two major parts. The first part contains alternating convolutional and pooling layers, which is used for feature learning. The second part is the fully connected layer, which can act as a classifier. The convolutional layer performs a dot product of the convolution kernel with the layer's input matrix. The pooling layer serves to progressively reduce the spatial size of the representation, to reduce the number of parameters and amount of computation in the
network, and hence to also reduce computation cost. There are several nonlinear functions to implement pooling, where max-pooling is the most common. It calculates the maximum value for patches of a feature map, and creates a downsampled (pooled) feature map. Meanwhile, this also implies that pooling operation would result in unacceptable information loss.

Recently, many optical neural networks based on various architectures have been proposed and experimentally demonstrated with impressive performance. However, these hardwares only support small matrices and implement few neurons. To make large neural networks run on the implemented photonic hardware accelerators, the down-sampling operations are usually used in the preprocess before photonic computing[39], and simplified neural network models are also utilized[42, 27, 37].

Figure 3: **An OCNN for image classification.****a** Schematic of OCNN-based image classification. The OCNN consists of a convolutional layer and a fully connected layer. The linear operations in these two layers are implemented by the photonic MAC operation unit. **b** The output temporal waveform for the feature map extracting. Inset shows the reconstructed feature map. **c** The output temporal waveform from the fully connected layer. The sampling points from left to right correspond to the handwritten digits from 0 to 9. The maximum sampling point corresponds to the classification result.
Here, based on the proposed photonic MAC operation unit, we further construct an OCNN for classifying handwritten digit images in the MNIST dataset. Thanks to its intrinsic advantage of temporal large-scale matrix operation, the max-pooling layer is not used in the designed CNN architecture. As shown in Fig. 3(a), a 5\(\times\)5 convolution kernel is used to extract the feature map of the 28\(\times\)28 handwritten digit images in the convolutional layer. By sampling the output temporal waveform of the photonic MAC operation unit and using the activation function (ReLu) to process the sampling result, a vector containing 784 elements is obtained. Then arranging these 784 elements into a 28\(\times\)28 matrix, the feature map of the input handwritten digit image is obtained in Fig. 3(b). The experimentally reconstructed feature maps from the convolutional layer exhibit good agreement with the theoretically calculated results (See Supplementary Fig. 5 and 6 for more detailed information). The photonic MAC operation unit can be further used to implement a fully connected layer, because of the linear matrix operation of the fully connected layer. To realize a complete functionality of an OCNN, the obtained feature maps from the convolutional layer are loaded into the photonic MAC operation unit again, which functions as a fully connected layer this time. In order to match the output of the convolutional layer, the fully connected layer is formed with 784 inputs and 10 outputs. Thus, the weight bank of the fully connected layer is a 10\(\times\)784 matrix. The weight bank matrix could be divided into 10 vectors with a length of 784, and sequentially loaded on the photonic MAC operation unit together with the feature map vectors. The output temporal waveform from the fully connected layer is recorded in Fig. 3(c). The sampling points from left to right correspond to the handwritten digits from 0 to 9. The image classification result is determined by the serial number of the maximum value in the output waveform of the fully connected layer. As can be seen in Fig. 3(c), the maximum sampling point corresponds to the handwritten digit 2, which agrees with the theoretically calculated one. We use the trained 5\(\times\)5 convolution kernel to process 20 handwritten images to get their feature maps, and the corresponding classification result gives a classification accuracy of 85% (see Supplementary Note 3 for more detailed information). In addition, more classification tests are performed by using 100 theoretically calculated feature maps. The confusion matrix for 100 images shown in Fig. 4(a) shows a classification accuracy of 88%, in contrast to the theoretically calculated accuracy of 93.2% on a digital computer (Fig. 4(b)). Due to the fluctuations of the device performances, the experimental accuracy will be slightly lower than the theoretically calculated value under the adopted CNN model. By integrating the discrete devices on a small single chip, the temporal computing system would be more stable.

In the experiment, both of the image data and the convolution kernel are flattened into temporal data groups. Each data group contains 25 symbols with a symbol rate of 6.92 gigabaud per second, and the time interval between each data group is 10 ns. The bandwidth of the BPD used for accumulation is 150 MHz. Therefore, in the convolutional layer, the total computing speed is 2\(\times\)25/10=5 GMAC operations per second (see Supplementary Note 6). In the implementation of the fully connected layer, the modulated symbol rate
is set to be 20.56 gigabaud per second, and a BPD with a bandwidth of 4 MHz is used to achieve the accumulation operation. As shown in Fig. 2c, the time interval between data groups is 583.66 ns, which could be further shortened. The effective computing speed for the fully connected layer is 2\(\times\)784/583.66=2.69 GMAC operations per second (see Supplementary Note 6).

As mentioned above, the down-sampling process and model simplification will cause some input information to be lost, thereby reducing the computational accuracy of neural networks. To investigate the impact of the down-sampling operation on the computational accuracy of neural networks, a simulated model shown in Fig. 3c is built in TensorFlow, which is trained by the AdamOptimizer with a default learning rate of 0.001, a training period of 50 iterations and a batch size of 1024. Same as the architecture used in the OCNN experiment, the convolutional layer contains a 5\(\times\)5 kernel and a nonlinear activation function (ReLu). A max-pooling layer is employed after the convolutional layer in the simulated CNN model. The training accuracies of the CNN model with four different down-sampling rates (1\(\times\), 4\(\times\), 9\(\times\), and

Figure 4: **Performance analysis of the proposed OCNN.****a** Confusion matrix obtained from 100 repeated experiments. **b** Confusion matrix obtained by simulation with 10,000 samples. **c** A common CNN architecture is composed of a convolutional layer, a pooling layers and a fully connected layer. **d** Sceneries are investigated. The training curves for CNN with different down-sampling rates are shown. Notably, the accuracy that the CNN model can obtain decreases with the increase of the down-sampling rate in the pooling layer.
16\(\times\)) are simulated. As can be seen in Fig. 4(d), the simulated classification accuracy that the CNN model can achieve becomes lower as the subsampling rate increases. Compared with the simulated accuracy of 93.2% in the CNN model without the down-sampling operation, 4\(\times\), 9\(\times\), and 16\(\times\) down-sampling will lead to training accuracies of 90.6%, 85.9%, and 80.7%, respectively. Hence, the CNN without pooling layer significantly outperforms its down-sampling counterparts. Similiar performance analysis for the fully connected neural network is carried out in Supplementary Fig. 10, which further verifies the feasibility and advantages of the proposed approach. The temporal computing approach holds great potential to enhance the learning capability of optical neural networks and processing more complex computational tasks, thereby helping to alleviate the problems imposed by the current limitations of manufacturing large-scale photonic integrated chips.

Figure 5: **Object detection based on sliding-window method.****a** The principle of object detection based on sliding-window method. **b** Schematic of using photonic MAC operation unit to detect the handwritten digits in a 68\(\times\)68 image. **c, d, e** The detection results of the handwritten digits 0, 4, and 8, the order of the detected peaks in the waveforms implies the position of the object in the original image.
**Sliding-window object detection with photonic MAC operation unit.** The sliding-window method is a very classical approach in object detection, and it also forms the basis of many other advanced object detection algorithms[48, 49]. Here, the photonic MAC operation unit is further used to realize a primary sliding-window multi-target detection. As shown in Fig. 5(a), the procedure of the sliding-window method is to crop the image into a series of image patches and classify all image patches as objects or non-objects independently. Figure. 5(b) shows the detail of using the photonic MAC operation unit to detect the handwritten digits in a 68\(\times\)68 image. The size of the sliding-window is set to be 28\(\times\)28 to match the size of the handwritten digits in the image. By moving the window in steps of ten pixels, the 68\(\times\)68 image is divided into 25 sets of 28\(\times\)28 image patches. All the 25 image patches are sequentially fed into a linear classifier, which is the most computationally intensive part of the sliding-window method. In the simulation, the linear classifier with a 28\(\times\)28\(\times\)10 three-dimensional (3D) matrix is trained in TensorFlow through the MNIST dataset. According to the classification results, the linear classifier matrix is flattened into 10 vectors with a length of 784. The image patches and the trained vectors are loaded into the photonic MAC operation unit to imitate the process of the classifier sliding through the image.

By sampling the output temporal waveform of the photonic MAC operation unit, the decision values for all the image patches are obtained. A decision threshold is adopted to judge whether there is an object in the image patch. In the experiment, the symbol rate is 20.56 gigabaud per second, and a time interval of 243.19 ns is used between each data group. The effective computing speed is 2\(\times\)784/243.19=6.45 GMAC operations per second. Figure. 5(c)-(e) show the corresponding output waveforms from the BPD for the handwritten digits 0, 4, and 8 shown in Fig. 5(b). The order of 25 pulses in the waveform represents the order in which the image patches are divided from the original image. The sample points that exceed the decision threshold determine the position of the object in the original image. For example, the output waveform of the handwritten digit 8 is shown in Fig. 5(e), where the value of the 5th sampling point from left to right exceeds the decision threshold. Thus, the object for the handwritten digit 8 is detected in the first row and fifth column of the original image. The object detection for other 4 sets of handwritten digits has also been successfully achieved in the experiment (see Supplementary Note 4). Here, a linear classifier is used instead of a multi-layer neural network classifier for multi-target detection detection. The straightforward approach adopted in the proposed model enables more than 4 million times objected or non-objected judgments per second in the designed task. The sliding-window object detection model could be further expanded to a more powerful model, such as a multi-layer neural network for more complex tasks such as face detection[50, 51].
Figure 6: **Parallel fully connected neural network for image recognition task.****a** Schematic of the parallel fully connected nueral network. The input image signals are multiplexed on two different light wavelength, and transmit through the implemented optical neural network. **b** The transmission spectrum of the multiplexed optical signals. **c,** The transmission spectrum of the multiplexed optical signal with one channel filtered out by the optical filter. **d, e** The measured waveforms for the handwritten digit 4 at 1549.8 nm and the handwritten digit 2 at 1550.2 nm, respectively. The maximum sampling point corresponds to the classification result. **f** The normalized sampling results at 1549.8 nm. **g** The normalized sampling results at 1550.2 nm.
**WDM-compatible parallel photonic computing.** In addition to the ultrahigh bandwidth inherent in optics, the ability of high parallel processing is another major advantage. Thanks to the scalability and the broadband operating characteristic of the proposed photonic MAC operation unit, it can be further expanded and enhanced by multiplexing wavelength and space dimensions. As a proof-of-concept experimental demonstration, we construct a WDM-compatible fully connected neural network, where two wavelength-multiplexed channels are employed to classify the handwritten digit images in the MNIST dataset. The experimental setup of the WDM-compatible fully connected neural network is shown in supplementary Fig. 9. Two different handwritten digital images with 28\(\times\)28 pixels are separately flattened into two vectors and loaded on two different wavelengths through two MZMs in the first stage. Then the temporal signals carried by two different light wavelengths are combined through an optical coupler as the input of the fully connected layer. The diagram of a parallel fully connected neural network is shown in Fig. 6. Two 1\(\times\)784 vectors flattened by two 28\(\times\)28 images are loaded on the first-stage MZM, and a 10\(\times\)784 matrix as the weight bank is loaded onto the second-stage MZM. After completing the signal multiplication operation, optical bandpass filters are used to filter out one of the wavelengths in Fig. 6(b) and (c). Then the filtered light signal is sent to the BPD for the temporal accumulation. The final step is to sample the temporal waveform to get the output results of the fully connected layer. As can be seen in Fig. 6(d) and (e), the maximum sampling points in the measured temporal waveforms correspond to the handwritten digit 4 at 1549.8 nm and the handwritten digit 2 at 1550.2 nm, respectively (see Supplementary Note 5 for more experimental results). The normalized sampling results for handwritten digits from 0 to 9 at the wavelength of 1549.8 nm and 1550.2 nm are summarized in Fig. 6(f) and (g), which agree well with the theoretically calculated ones. In the experiment, the symbol rate is 20.56 gigabaud per second, and a time interval of 583.66 ns is used between each data group.Thus, the effective computing speed is 2\(\times\)784/583.66=2.69 GMAC operations per second for one wavelength. The total computing speed of the WDM-compatible fully connected neural network is 2\(\times\)2.69=5.38 GMAC operations per second (see Supplementary Note 6). The proposed computing architecture can significantly enhance the computation capability of optical neural networks with more multiplexed channels.
## Discussion

We propose a compact and scalable photonic computing architecture based on temporal photoelectric MAC operations, where large-scale matrix multiplication operations could be directly processed in the time domain. Particularly, the proposed approach offers the opportunity to extend the length of operand vectors and avoid additional down-sampling process in neural networks, which will be helpful to obtain higher computing accuracy in a range of tasks. Based on the photonic MAC operation unit, a series of proof of principle experiments are demonstrated, including high-speed convolutional image processing, optical neural networks-based recognition tasks, and sliding-window method-based multi-target detection. Furthermore, based on the proposed parallel photonic computing architecture, we construct a WDM-compatible fully connected neural network for image recognition task.

Recently, a lot of approaches have been experimentally demonstrated to realize photonic matrix multiplication and nueral networks. Table 1 shows the comparison between the reported works and our work, and some key metrics evaluating the physical network are listed, including the architecture, accumulation mechanism, and the matrix dimenison. Previously MZI-based integrated photonic approaches have been predominantly limited by the current limitations of manufacturing large-scale photonic integration hardwares, resulting in large footprints and small matrix sizes (4\(\times\)4 and 6\(\times\)6). Addtionally, the accumulation operation of MRR/PCM array-based solutions is accomplished by summing up the light fields

\begin{table}
\begin{tabular}{c c c c c} \hline
**Year** & **Ref.** & **Accumulation approach** & **Architecture** & **Matrix dimension** \\ \hline
2017 & 25 & Coherent light accumulation & MZI array & 4\(\times\)4 \\
2021 & 27 & Coherent light accumulation & MZI array & 6\(\times\)6 \\
2021 & 28 & Coherent light accumulation & MZM array & 9\(\times\)1 \\
2019 & 33 & Multi-wavelength accumulation & MRR \& PCM array & 4\(\times\)15 \\
2021 & 36 & Multi-wavelength accumulation & MRR array & 4\(\times\)2 \\
2021 & 37 & Multi-wavelength accumulation & PCM array & 16\(\times\)16 \\
2021 & 40 & Multi-wavelength accumulation & Discrete device & 90\(\times\)1 \\
2021 & 43 & Multiple photodetectors accumulation & MRR \& VOA array & 5\(\times\)6 \\
2022 & 53 & Coherent light accumulation & Cascaded MMI & 10\(\times\)10 \\  & Ours & Time-domain integration & MZM array & 784\(\times\)1 \\ \hline \end{tabular}
\end{table}
Table 1: Comparison between the reported optical neural networks and our work.
of multi-wavelength channels into a single photodetector, which also leads to a limited matrix dimension (16\(\times\)16). And other approaches are also limited by either the integration scale or the number of wavelength channels. By comparison, our approach is based on time-domain integration, which determines that the length of the vector for dot product is therectically unlimited. In the experiment, the length of the loaded vector can be up to 784, which is much larger than that of previous approaches. Furthermore, the demonstrated photonic computing unit could be easily integrated on a chip and scaled up to a multidimensional parallel photonic MAC operation architecture by employing time, wavelength and spatial dimensions, which can further increase the number of parallel computation channels. The advantages are briefly couccluded as: (a) It offers the ability to perform large-scale matrix multiplication operation, which is determined by the temporal MAC operation mechanism. (b) Its strong scalability and high parallelism make it possible to use multiple dimensional resources simultaneously for computational acceleration.

In this work, the effective computational speed is performed at the level of Giga MAC operations per second, which is restricted by the modulators and BPD. With further improvement, the state-of-the-art integrated photonic transmitters and photodetectors can drive the MAC operation unit at a speed of 100 GOPS[52]. When scaling up to a larger network with 100 parallel wavelength/space computing channels, the total effective computing speed can easily reach 10 Tera MAC operations per second (TOPS). Thus, the proposed photonic computing architecture is capable of performing large-scale matrix-vector multiplication. Furthermore, it is worth mentioning that the proposed temporal computing architecture can be easily fabricated on a small single chip, and thus its enhanced functionality and performance will make it practical to implement large-scale matrix computation for photonic neural networks. Additionally, complex-valued neural network could be realized by IQ modulation and receiver, which would further improve the computing performance.

## Methods

### Experimental setup.

The light source used in the experiments was NKT Koheras BASIK 1550 nm single-frequency fiber laser. A 4-channel Keysight M8196A AWG with 92 GSa/s, 8-bit vertical resolution and 32 GHz analog bandwidth was used for digital-to-analog converting. The electro-optic modulation was all conducted by Fujitsu FTM7937 modulators. The BPD is Thorlabs PDB 450C with switchable gain and bandwidth. A Tektronix DPO75002SX Oscilloscope was used to analyze the output waveform of BPD. The digital filter of the oscilloscope was set to 500 MHz to suppress the high-frequency noise of the BPD.