DT-OCS is conceived as a system-level digital twin framework that models the end-to-end input-output behavior of a physical OCS. To illustrate the principles and construction methodology of DT-OCS, we instantiate the framework on a representative high-speed silicon-based integrated optical computing system, shown in Fig. 2(a), as the physical platform for validation. The optical computing chip in this system has high-speed temporal signal processing capabilities and supports programmable computing functions. Our previous work has validated that this type of physical OCS can compute the feature information of high-speed signals^33,34. As shown in Fig. 2(c–e), the optical computing chip within the system consists of two core modules. 1) Data preprocessing module: This module performs sampling, serial-to-parallel conversion, and phase adjustment on the time-sequenced input signal $ {{{{{i}}}}}(t) $, generating spatially parallel signalwhere $ {\mathrm{e}}^{{{\mathrm{i}}}\varphi_n} $ and $ {{{\Delta}} }_{n} $ represent the added phase and delay time, respectively. By designing the lengths of the delay waveguides, the sampling rate is set to 10 GHz. 2) Feature extraction core: This core performs matrix multiplication. Based on the input vector $ {\boldsymbol{I}}(t) $, the feature extraction core computes the output $ O(t)={{\boldsymbol{SI}}}(t) $, where $ {\boldsymbol{S}\mathit{ }} $ represents the port scattering matrix of the feature extraction core. Under this operating condition, the physical OCS performs computation directly at the symbol rate, such that the computation itself does not introduce additional clocked digital latency. Using this optical computing chip package, the estimated throughput^25,28,35 and equivalent energy efficiency³⁶ reach 100 GOPS and 0.8 TOPS/W, respectively. Details of the optical computing chip initialization and the estimation of system-level metrics of the physical OCS are provided in Supplementary Section 1 and Section 2, respectively.

However, the actual output deviates from the expected behavior: firstly, parasitic dispersion, reflection, and resonance in the feature extraction core and other components in the waveguide optical path introduce wavelength-dependent responses. As a result, the output of the physical OCS is influenced by an extended input receptive field, which not only samples more discrete points but also incorporates adjacent time ranges. Similarly, as shown in Fig. 2(a), signal-related devices in the physical OCS further complicate this issue. These devices include the signal source device (AWG), modulation amplifier, phase modulator (PM), optical amplifier (EDFA), PD, and sampling device (OSC), all of which have limited bandwidth or non-ideal responses. Together, these factors result in a complex input-output mapping of the physical OCS that deviates from what can be described by simple linear matrix multiplication, necessitating a more generalized and systematic modeling approach to study the physical OCS. Detailed information on the system’s components can be found in Supplementary Section 2.

The physical OCS is a finite memory system, whose output at a given time $ t $ depends on the input from a past time period, i.e.,

The operator $ \mathcal{F} $ represents the computational functionality of the physical OCS and captures the overall system response resulting from the optical processor and the associated components, thereby defining the mapping rule between input and output data. The input $ {{i}}(\xi ) $ over a time period $ [t-\tau ,t] $ determines the system's output, where $ \tau $ represents the time length of the receptive field. $ \boldsymbol{{V}} $ represents the bias voltage of the phase shifters in the data preparation module of the optical computing chip in the physical OCS. By changing the bias voltage (configuration parameter) $ {\boldsymbol{V}}=[V_p] $, the additional phase in the branch is altered, thus adjusting the computational functionality of the physical OCS. Here $ p $ represents the number of voltage configuration parameters, which corresponds to the number of actively controlled phase shifters in the OCS. The computational functionality $ \mathcal{F} $ of the physical OCS is realized only through direct hardware involvement, and becomes inaccessible when the physical system is offline. By constructing DT-OCS, we create an open-source and callable operator $ \hat{\mathcal{F}} $ as a surrogate for the physical operator $ \mathcal{F} $, enabling DT-OCS to participate in task design and training while the physical OCS is offline.

The input signal to the physical OCS is a piecewise constant signal generated from the data to be computed. Therefore, the vector $ {\boldsymbol{I}}_{{{t}_{0}}}=[{{i}}({t}_{0}-(K-1)\Delta ), \cdots ,{{i}}({t}_{0}-k\Delta ), \cdots , {{i}}({t}_{0})] $ is able to encapsulate the input information for the physical OCS within the receptive field $ [{t}_{0}-\tau ,{t}_{0}] $. The number of elements in vector $ {\boldsymbol{I}}_{{{t}_{0}}} $ is denoted aswhere $ \mathrm{ceil}[\cdot ] $ is the ceiling function, $ {{\Delta}} $ represents the modulation symbol duration, $ \tau $ denotes the receptive field length, and $ {t}_{0} $ refers to the sampling completion time point. The response of the physical OCS is similarly sampled at the same rate to form $ {O}_{{{t}_{0}}} $. In this way, the computational function of the physical OCS can be described as follows: under configuration parameter $ \boldsymbol{V} $, the physical OCS generates an output $ {O}_{{{t}_{0}}}=\mathcal{F}({\boldsymbol{I}}_{{{t}_{0}}};\boldsymbol{V}) $ based on the input $ {\boldsymbol{I}}_{{{t}_{0}}} $.

According to the finite memory and causal characteristics of the physical OCS, we adopt a long short-term memory (LSTM) architecture³⁷ as the digital foundation of DT-OCS. To embed the configuration-dependent computational function of the physical OCS, the standard LSTM is modified such that its computational behavior is governed by the same configuration parameters $ \boldsymbol{V} $. Specifically, as shown in Fig. 2(b), the internal parameters, including the weight matrices and bias vectors associated with the LSTM gates, are not fixed but are generated by a function block that takes the physical configuration parameters $ \boldsymbol{V} $ as input. This function block is realized as a simple linear (affine) mapping from $ \boldsymbol{V} $ to the gate parameters of the LSTM, while all matrices and bias in this mapping are trainable and collectively denoted as $ \theta $. In this manner, a standard LSTM cell is transformed into a Dynamic LSTM cell, whose recurrent dynamics explicitly depend on the configuration parameters of the OCS. The number of cells in DT-OCS is also set to $ K $, ensuring that DT-OCS shares the same receptive field as the physical OCS. The final output of DT-OCS is $ {\hat{O}}_{{{t}_{0}}}=\hat{\mathcal{F}}({\boldsymbol{I}}_{{{t}_{0}}};\boldsymbol{V},\theta ) $, where $ \hat{\mathcal{F}} $ is realized with the DT-OCS network. DT-OCS is expected to be equivalent to the physical OCS, meaning that $ \hat{\mathcal{F}} $ closely approximates $ \mathcal{F} $. This equivalence is reflected in the fact that, given the same input $ \boldsymbol{I} $ and configuration parameters $ {\boldsymbol{V}} $, the output of DT-OCS $ \hat{O} $ is sufficiently close to the output of the physical OCS $ O $. This closeness is quantified by a similarity function (SF), which evaluates the mean squared error (MSE) between the two output sequences $ \hat{O} $ and $ O $ in both time and frequency domains:

Joint evaluation of errors in both time and frequency domains helps avoid degenerate similarity evaluations. Such degeneracy may arise when relying on a single representation, which tends to emphasize either local waveform details or global structural characteristics alone. Similar joint loss functions have been widely adopted in signal processing research^38,39. Detailed implementation of DT-OCS is available in the open-source code linked at the end of the paper.

We constructed a learning dataset for DT-OCS under the supervised learning framework. The learning objective of DT-OCS is to accurately approximate the input-output mapping of the physical OCS across different configuration parameters, as shown in Fig. 3(a). The dataset is organized in a triplet structure, where each training sample consists of a voltage configuration $ {\boldsymbol{V}}^{m} $, an input segment $ {\boldsymbol{I}}_{\{{{t}_{n}}\}} $, and the measured physical OCS output $ O_{\{{t}_{n}\}}^{m} $ under that configuration parameter, i.e., mapping pairs $ {\boldsymbol{V}}^{m} $−$ {\boldsymbol{I}}_{\{{{t}_{n}}\}} $−$ O_{\{{t}_{n}\}}^{m} $. Specifically, 4000 distinct voltage configurationsare loaded onto the phase shifters shown in Fig. 2(c) of the optical computing chip in the physical OCS, where $ m $ denotes the configuration parameter index. These configurations $ {\boldsymbol{V}}^{m} $ randomly sampled from a range spanning twice the half-wave voltage of the phase shifters. Under each configuration parameter, a fixed probing sequence is input into the physical OCS. This sequence consists of multimodal data, including image information, time series, and random noise, with a total length of 2500 points. Details of this sequence construction are provided in Supplementary Section 4. For a given configuration parameter, the physical OCS produces a corresponding output sequence in response to the input sequence. Representative segments of the input and the corresponding output sequence are shown in Fig. 3(f) and the blue markers in Fig. 3(g), respectively. Based on the each input-output sequence pair under each configuration parameter, a sliding time window is applied to the input sequence to extract input segments $ {\boldsymbol{I}}_{{{t}_{n}}} $ of length $ K $ ($ K $ is the element number of the input vector, as defined in Principles and methods). Each input segment $ {\boldsymbol{I}}_{{{t}_{n}}} $ is combined with the current configuration $ {\boldsymbol{V}}^{m} $ and the corresponding output $ O_{{t}_{n}}^{m} $ to form a complete mapping pair. Here, $ {t}_{n}\in \{K,K+1, \cdots, 2500\} $ denotes the sampling time. To determine the system’s receptive field $ \tau $ and accordingly compute $ K={\mathrm{ceil}}[\tau /\Delta ] $, a full-amplitude step function, shown in Fig. 3(b), is used to excite the physical OCS. The decay of the output response is recorded in Fig. 3(c). We define the receptive field $ \tau $ as the time interval over which the effect output response (the DC component 0.5 subtracted) decays from its maximum to 15%, which is approximately 0.73 ns. The 15% level is adopted as an operational criterion corresponding to a point at which the output response has essentially converged, leaving only a weak tail contribution. Given the system's 10 GHz sampling rate ($ {{\Delta}} =0.1\;\mathrm{ns} $), we determined $ K=8 $. Under this configuration, each voltage setting yields $ 2500- K=2492 $ mapping pairs, resulting in a total of 2492 × 4000 $= $ 9.968M mapping pairs, which form the supervised learning dataset for DT-OCS: $ {\boldsymbol{V}}^{m} $−$ {\boldsymbol{I}}_{\{{{t}_{n}}\}} $−$ O_{\{{t}_{n}\}}^{m} $. The range of the half-wave voltage is provided in Supplementary Section 1.

The learning of DT-OCS is then conducted. The learning objective is for the computational function of DT-OCS $ \hat{\mathcal{F}}(\boldsymbol{I};\boldsymbol{V},\theta ) $ to accurately approximate the function of the physical OCS $ \mathcal{F}(\boldsymbol{I};\boldsymbol{V}) $. The previously proposed time-frequency composite similarity function ($ \mathrm{SF} $) is selected as the loss function and optimization criterion for training the DT-OCS. By minimizing $ \mathrm{SF} $ via gradient descent, DT-OCS is guided to progressively approximate the input-output mapping of the physical OCS across the entire configuration parameter space. During the learning process, the configuration parameters $ \boldsymbol{V} $ are first passed into the function block defined by parameters $ \theta $, which then generates the internal parameters of the Dynamic LSTM cells. The learning process thus achieves the optimization of the function block parameters:

The 4000 voltage configurations are randomly divided into a train set and a test set in a 4∶1 ratio. As shown in Fig. 3(d), after 200 epochs of training, $ {SF} $ converges sufficiently. From the test set, we select one configuration and compare the outputs of DT-OCS and the physical OCS, as shown in Fig. 3(g). DT-OCS accurately fits the output of the physical OCS. Across the entire test set, the residuals between $ {\hat{O}}_{n} $ and $ {O}_{n} $ are fitted using a Voigt profile. The standard deviation of its Gaussian component is $ {\sigma }_{\mathrm{Gaussian}}=3.24\times {10}^{-3} $, and that of its Lorentzian component is $ {\sigma }_{\mathrm{Lorentzian}}=\mathrm{4.82}\times {10}^{-3} $. Both are below 0.5‰ of the maximum output amplitude, demonstrating the accuracy of the fit. At this point, the learning process of DT-OCS is complete, and the trained DT-OCS is used for subsequent training and testing tasks.

The residual between the output of the physical OCS and that of its digital twin model, DT-OCS, serves as a key metric for evaluating both the computational reliability of the OCS and the accuracy of the DT-OCS framework. The origins of the residual can be attributed to two primary sources: 1) Model approximation error. DT-OCS approximates the computational function of the physical OCS by mapping voltage configuration $ \boldsymbol{V} $ to the internal parameters of the Dynamic LSTM through a function block. In the current implementation, this function block adopts a single-layer linear mapping, as illustrated in the inset of Fig. 7(a). The resulting residual distribution on the test set, together with its Voigt profile fitting shown by the red curve, is shown in Fig. 7(a) (identical to Fig. 3(e)). While this linear design offers high computational efficiency and stable training behavior, it may not fully capture nonlinear coupling mechanisms among configuration parameters in the physical OCS. A representative example is thermal crosstalk among thermo-optic phase shifters, where localized temperature variations induced by different voltage configurations lead to heat diffusion and mutual interference, giving rise to nonlinear parameter coupling that is difficult to represent using a purely linear mapping^41,42. To examine the impact of enhanced parameter mapping, we introduce a more complex function block by inserting an intermediate nonlinear transformation between two linear layers, as illustrated in Fig. 7(b). The detailed layer architecture and parameter dimensions of both the linear and nonlinear function blocks are explicitly provided in the Supplementary Section 3. Using this enhanced function block, we re-evaluate the residual distribution and again fit the distribution with a Voigt profile. The standard deviation of the Gaussian component is $ {\sigma }_{\mathrm{Gaussian}}=3.41\times {10}^{-3} $, and that of the Lorentzian component is $ {\sigma }_{\mathrm{Lorentzian}}=\mathrm{3.25}\times {10}^{-3} $. Compared to the single-layer configuration, the standard deviation of Gaussian component of the fitted distribution increases by 5.24%, while the standard deviation of Lorentzian component decreases by 32.57%, indicating a more Gaussian-like distribution with suppressed long-tail behavior. As shown by the yellow curve in Fig. 7(b), the residuals obtained with the nonlinear complex function block are more concentrated than those of the simple linear design (red dashed line), demonstrating improved fitting fidelity. However, this improvement comes at the cost of increased computational complexity. A quantitative comparison of the associated parameter counts and computational complexity is provided in the Supplementary Section 3. To balance modeling accuracy and efficiency, the simplified linear function block is adopted as the primary configuration in this work, while the nonlinear alternative is also provided for readers interested in exploring this trade-off further. A quantitative comparison of fitting statistics for both designs on the training and test sets is reported in Supplementary Section 3. Implementations of both designs are provided in the open-sourced codebase referenced in the Code Availability section. 2) Physical noise from hardware. Residuals also arise from stochastic noise sources inherent to the physical OCS hardware. Signal generation, modulation, and readout components, including the AWG, modulation amplifiers, high-speed phase modulators, PD, and associated electronic circuits, introduce timing jitter, thermal noise, and random fluctuations. Instabilities in waveguides and coupling fibers can further lead to phase noise, while optical amplification contributes to gain noise. These physical noise sources often manifest as statistically distributed disturbances, yet they are difficult to capture accurately using a single general model or fixed functional form. In summary, the residuals arise from a combination of model approximation error under practical design choices and hardware-induced physical noise mechanisms. Figure 7(a) shows that the Voigt profile (red line) matches the residual distribution more closely than either a purely Gaussian (green line) or Lorentzian (blue line) profile. This likely results from a combination of the aforementioned factors. Although the architecture and operation differ, MZI mesh-based structures in coherent optical systems have exhibited similar residual behavior⁶. Once the residual distribution is obtained, noise following the same distribution can be injected into the DT-OCS output via rejection sampling. This offers a flexible mechanism to incorporate measured noise into the digital twin framework when robustness evaluation or data augmentation is desired. The concrete implementation details of the noise-generation module are provided in the Code Availability section.

To investigate potential nonlinear response characteristics that may arise during OCS operation, we design an excitation experiment with inputs of identical shape but varying amplitude, as illustrated in Fig. 8(a). The system response is analyzed using DT-OCS. Specifically, as shown in Fig. 8(b), we fix the shape of the input signal, a threshold signal with a width of 8, and gradually increase its amplitude by adjusting the amplification factor Amp, over the range from 0 to 1. The configuration parameter is selected from the supervised learning dataset ($ {\boldsymbol{V}}^{2000} $).

In an ideal linear system, responses to inputs with the same shape but different amplitudes should completely overlap after bias removal and normalization (i.e., division by the amplification factor Amp). However, as shown in Fig. 8(c), after the same bias removal and normalization, the DT-OCS output exhibits significant variations not only in amplitude but also in shape and peak position as the input amplitude increases. Additionally, the gray lines in Fig. 8(c) indicate the DT-OCS outputs at two specific positions P=10 and P=18, whose raw bias-removed outputs as functions of the input amplitude are shown by the colored curves in Fig. 8(d) and 8(e), respectively. As the input amplitude increases, the two positions show markedly different nonlinear behaviors in their raw output amplitudes. The same experiment is repeated on the physical OCS, where the corresponding raw bias-removed outputs are shown as discrete points in Fig. 8(d) and 8(e). Despite measurement noise and deviations, the physical OCS outputs exhibit nonlinear response trends consistent with those of the DT-OCS, with the two positions again showing markedly different nonlinear behaviors. These observations collectively demonstrate the presence of complex internal nonlinear mechanisms within the system. Possible sources of the observed nonlinearity include phase-to-intensity mapping, non-ideal behavior of high-speed electro-optic modulators, and saturation effects in the modulation amplifier. DT-OCS captures and reproduces these nonlinear response characteristics, enabling system-level modeling for task training and parameter optimization.

Accurate modeling and deployment of high-speed OCSs require a well-defined timing relationship between the input and output signals, which in practice is established through system-level synchronization between signal generation and acquisition. In our experimental setup, the AWG serves as both the signal source and the timing reference. It generates the input signal and provides a trigger signal to the oscilloscope, establishing stable synchronization during signal acquisition. To provide a consistent temporal reference across measurements, a predefined threshold-based header pattern is inserted before each data sequence and used as the trigger signal for the oscilloscope. The rising edge of this header acts as a temporal alignment marker. In the recorded output signal, the corresponding initial rise of the measured output response indicates the onset of the system's causal response, consistent with the beginning of the receptive field marked by the left gray dashed line in Fig. 3(c), thus establishing a well-defined input-output correspondence. By fixing this reference point, the portion of the output response corresponding to the data payload in the header-augmented input signal can be consistently identified, enabling reliable subsequent segmentation of the continuous output signal into valid output sequences. An illustrative example is provided in Supplementary Section 6. This input-output correspondence is determined during the receptive-field measurement stage and is subsequently maintained throughout data acquisition, training, and inference, ensuring causal consistency between the physical OCS and its digital twin. This mechanism effectively avoids causal ambiguity in high-speed OCS operation and forms the basis for reliable DT-OCS modeling.

In this work, DT-OCS is instantiated on a representative high-speed silicon-based optical computing system to validate the feasibility of the proposed framework. This implementation demonstrates that the DT-OCS methodology can be constructed and deployed on a practical optical computing platform. DT-OCS models the end-to-end behavior of the assembled system under a given hardware configuration, thereby reflecting the system-level characteristics of the underlying physical platform.

In practical operation, hardware configuration may change due to component replacement, system upgrade, or architectural modification. The impact of such changes on DT-OCS depends on the scale and structural level of the modification. When optoelectronic components are replaced by devices with the same specifications and operating characteristics, the system-level input-output behavior is expected to remain consistent, and the trained DT-OCS remains applicable. If replacement introduces changes in device-level characteristic parameters, such as bandwidth, noise level, gain, or nonlinear response, that in turn alter the system-level input-output behavior, the digital twin can be updated by re-measuring system-level input-output data and re-training the DT-OCS while retaining the same modeling framework.

When the optical processor is redesigned into a different computational architecture, the underlying computational function changes, and the previously trained DT-OCS instance no longer corresponds to the updated system. In this case, as long as the updated optical computing system defines a causal input-output mapping governed by configuration parameters, a corresponding DT-OCS can be constructed under the same measurement-driven framework. Although different platforms may require adjustments in modeling details or training procedures due to variations in physical implementation, the overall methodology and workflow remain unchanged. In this sense, DT-OCS should be understood as a reproducible framework for encapsulating the configuration-dependent behavior of an optical computing system into a digital operator, supporting systematic modeling and integration of optical processors within hybrid optoelectronic computing systems.