Role of Spatial Coherence in Diffractive Optical Neural Networks
Diffractive optical neural networks (DONNs) have emerged as a promising
optical hardware platform for ultra-fast and energy-efficient signal processing
for machine learning tasks, particularly in computer vision. Previous
experimental demonstrations of DONNs have only been performed using coherent
light. However, many real-world DONN applications require consideration of the
spatial coherence properties of the optical signals. Here, we study the role of
spatial coherence in DONN operation and performance. We propose a numerical
approach to efficiently simulate DONNs under incoherent and partially coherent
input illumination and discuss the corresponding computational complexity. As a
demonstration, we train and evaluate simulated DONNs on the MNIST dataset of
handwritten digits to process light with varying spatial coherence.
Matthew J. Filipovich, Aleksei Malyshev, A. I. Lvovsky
# Role of Spatial Coherence in Diffractive Optical Neural Networks

###### Abstract

Diffractive optical neural networks (DONNs) have emerged as a promising optical hardware platform for ultra-fast and energy-efficient signal processing for machine learning tasks, particularly in computer vision. However, previous experimental demonstrations of DONNs have only been performed using coherent light, which is not present in the natural world. Here, we study the role of spatial optical coherence in DONN operation. We propose a numerical approach to efficiently simulate DONNs under input illumination with arbitrary spatial coherence and discuss the corresponding computational complexity using coherent, partially coherent, and incoherent light. We also investigate the expressive power of DONNs and examine how coherence affects their performance. In particular, we show that under fully incoherent illumination, the DONN performance cannot surpass that of a linear model. As a demonstration, we train and evaluate simulated DONNs on the MNIST dataset of handwritten digits using light with varying spatial coherence.

## I Introduction

Machine learning has led to breakthroughs in many fields of technology including computer vision, natural language processing, and drug discovery [1, 2, 3, 4]. Currently, machine learning models are executed using specialized electronic hardware, such as graphics processing units (GPUs) and tensor processing units (TPUs), which harness immense processing power and data parallelism. However, the growing compute requirements from advanced deep learning models are far outpacing hardware improvements anticipated by Moore's law scaling [5]. Consequently, progress in machine learning using current digital hardware will soon become technically and economically unsustainable, while also producing a significant environmental impact [6, 7].

Given the constraints imposed by digital electronics, optics has gained recognition as a promising platform for machine learning applications with low latency, high bandwidth, and low energy consumption [8]. Several optical implementations of neural networks have recently been demonstrated using both free-space optics [9, 10, 11, 12] and integrated photonics [13, 14, 15, 16]. Computational speeds of trillions of operations per second have been achieved using optical processing hardware [17, 18], and optical neural networks using less than one photon per multiplication have been experimentally realized [19, 20]. Additionally, optical architectures for implementing _in situ_ training of neural networks have been demonstrated [21, 22, 23, 24, 25, 26].

Diffractive optical neural networks (DONNs) are specialized hardware architectures that harness diffraction effects to process optical signals in free space [27, 28]. DONNs are generally composed of several successive modulation surfaces, denoted as diffractive layers, that modify the phase and/or amplitude of the incident optical signals through light-matter interactions, as shown in Fig. 1a. The diffractive layers contain discrete pixels, each with an independent complex-valued transmittance coefficient. The output of the DONN corresponds to the total intensity of the optical field incident on designated detection regions in the output plane.

This architecture is particularly promising for "deep optics" applications [29], i.e., complex processing and recognition of images without converting the optical signal into electronic digital form. In addition to machine learning tasks, DONNs have further applications in super-resolution imaging and displays, microscopy, quantitative phase imaging, and image classification through random diffusers [30, 31, 32, 33]. Several physical realizations of diffractive layers have been experimentally demonstrated using 3D printed materials, metamaterials, and spatial light modulators [34, 35, 36].

DONNs are usually _in silico_, i.e., the physical DONN is modeled using a computer to simulate the evolution of the optical signals through the system. The modulation patterns of the diffractive layers are optimized to achieve the desired transformation between the input and target output of the DONN, which is analogous to optimizing the weights in standard neural network models [37]. During training, the transmittance coefficient of each pixel in the diffractive layers is iteratively updated using an optimization algorithm to minimize the error in the model's output with respect to the training set. The backpropagation algorithm is used to efficiently calculate the gradient of the loss with respect to the transmittance coefficients [38].

DONNs are particularly well-suited for use with real-world optical signals, as the optical fields can be directly fed into the system. However, such signals are typically incoherent. In contrast, most of the existing experimental work with DONNs is performed using fully-coherent illumination from laser sources.

In this paper, we introduce a computationally efficient framework for simulating and training DONNs using incident illumination with arbitrary spatial coherence and discuss the corresponding computational complexity. We also investigate the role of optical coherence on the expressive ability of DONNs to extract complex relationships in the input data. In particular, we show that under fully incoherent illumination, the DONN performance cannot surpass that of a single-layer linear neural network. In contrast, when some degree of coherence is present, measuring the intensity at the output constitutes a nonlinear operation, which plays the role of an activation function in standard neural network models. This permits reaching performance
levels (e.g., classification accuracy) beyond the capabilities of a linear model. To illustrate our findings, we evaluate the performance of simulated DONNs trained on the MNIST dataset of handwritten digits using light with varying optical coherence.

A deep-learning based method was recently proposed for implementing linear transformations with DONNs under spatially incoherent illumination [39]. The method approximates incoherence for a single input example by averaging the output intensity patterns from numerous coherent input fields with random phase distributions. In contrast, our approach operates with the mutual intensity function (i.e., mutual coherence function), which is a compact and efficient way to represent the statistical properties of a partially coherent field [40]. This enables us to compute the propagation of light with arbitrary coherence through diffractive layers much more efficiently and accurately.

## II Results

### Coherent illumination

In this section, we introduce a formalism for describing the evolution of coherent, monochromatic optical fields through DONNs using scalar diffraction theory [41]. We treat the optical field as a complex scalar quantity and employ Dirac notation to represent the transverse profile of the field at discrete spatial positions using ket-vectors. This discretization does not affect the generality of our treatment as long as the spatial sampling interval (i.e., pixel pitch) is much smaller than the characteristic transverse field feature size. The transformations applied to the field as it evolves through the DONN, which include free-space propagation and transmission through modulation surfaces, are expressed using linear operators.

At each layer in the DONN, the optical field is modulated by the diffractive surface and subsequently propagates through free space to the next layer. The incident field at the \(m\)-th discrete pixel of the \(l\)-th diffractive layer, before modulation, is represented by \(\psi_{l}(m)\), where the time dependence of the signal is absent. The transverse profile of the field can be expressed using Dirac notation as \(\ket{\psi_{l}}=\sum_{m}\psi_{l}(m)\ket{m}\), where the set \(\{\ket{m}\}\) of all pixels forms an orthonormal basis. The mapping between the optical fields in the \(l\)-th and \((l\)+1\()\)-th layers can be expressed as \(\ket{\psi_{l+1}}=\hat{P}_{l}\hat{T}_{l}\ket{\psi_{l}}\), where \(\hat{T}_{l}\) and \(\hat{P}_{l}\) are the transmission and free-space propagation operators, respectively.

The transmission operator \(\hat{T}_{l}\) describes the phase and/or amplitude modulation applied to the optical field by each pixel in the \(l\)-th diffractive layer. The corresponding matrix is diagonal:

\[\hat{T}_{l}=\sum_{m}t_{l}(m)\cdot\ket{m}\!\!\bra{m}, \tag{1}\]

where \(t_{l}(m)\) is the complex-valued transmittance coefficient

Figure 1: Diffractive optical neural network (DONN) architecture. **a** Illustration of a DONN trained to identify handwritten digits. The DONN is comprised of \(L\) diffractive layers that modulate the optical field as it propagates through the system. The output plane encompasses ten detection regions, which are each associated with a unique digit, and the predicted output corresponds to the region with the highest optical intensity. The transmission and propagation operators at the \(l\)-th layer are denoted by \(\hat{T}_{l}\) and \(\hat{P}_{l}\), respectively. **b** Summary of the transformations applied to the input field \(\ket{\psi_{\mathrm{in}}}\) in a DONN. The linear evolution operator \(\hat{U}\) maps the input field onto the output field. The intensity \(I_{\mathrm{out}}\) at the output plane is measured, and the output of the DONN, represented by \(o\), corresponds to the total intensity incident on each detection region.
at the \(m\)-th pixel in the \(l\)-th diffractive layer, which satisfies \(|t_{l}(m)|\leq 1\).

The operator \(\hat{P}_{l}\) describes the free-space propagation of the field between the \(l\)-th and (\(l\)+1)-th diffractive layers using the Rayleigh-Sommerfeld solution [41]:

\[\hat{P}_{l}=\sum_{m,n}h(m,n)\cdot|n\rangle\!\langle m|\,, \tag{2}\]

with

\[h(m,n)=\frac{1}{i\lambda}\exp\left(\frac{i2\pi r(m,n)}{\lambda}\right)\frac{d }{r(m,n)^{2}} \tag{3}\]

being the point-spread function, i.e., the amplitude distribution in the (\(l\)+1)-th layer if only the \(m\)-th pixel of the \(l\)-th layer is illuminated. In the above equation, \(\lambda\) is the wavelength of the coherent optical signal (the central wavelength for quasimonochromatic light), \(d\) is the axial distance between the two diffractive layers, and \(r(m,n)\) is the Euclidean distance between the \(m\)-th and \(n\)-th pixels in the \(l\)-th and (\(l\)+1)-th layers, respectively. The above expression is valid when the axial distance between layers is much greater than the (central) wavelength of light.

The input image processed by the DONN is encoded in the initial field \(|\psi_{\text{in}}\rangle\). The output optical field, represented by \(|\psi_{\text{out}}\rangle\), can then be expressed as

\[|\psi_{\text{out}}\rangle=\hat{U}\,|\psi_{\text{in}}\rangle\,, \tag{4}\]

where \(\hat{U}\) is the evolution operator of the DONN that maps the input optical field onto the output field (i.e., the spatial impulse response of the system):

\[\hat{U}=\prod_{l=1}^{L}\left(\hat{P}_{l}\hat{T}_{l}\right)\hat{P}_{0}, \tag{5}\]

where the DONN has \(L\) diffractive layers. At the output plane of the DONN, the intensity of the evolved field is measured using image sensors:

\[I_{\text{out}}(n)=|\,\langle n|\psi_{\text{out}}\rangle\,|^{2}=|\psi_{\text{ out}}(n)|^{2}. \tag{6}\]

For classification tasks, the output of the DONN corresponding to each class \(c\) is defined as the total intensity incident on a specified spatial detection region \(\mathcal{D}_{c}\) in the output plane:

\[o(c)=\sum_{n\in\mathcal{D}_{c}}I_{\text{out}}(n). \tag{7}\]

A summary of the mathematical operations applied to the input optical signal by the DONN during inference is shown in Fig. 1b.

Training DONNs using a computer (i.e., _in silico_) requires simulating the evolution of coherent optical fields through the system. The calculated optical field at each layer is then used during the backward pass to compute the gradient of the loss function with respect to the diffractive layer transmittance coefficients. A naive numerical implementation of the propagation operator \(\hat{P}_{l}\) using matrix-vector multiplication has a computational complexity of \(\mathcal{O}(N^{2})\), where \(N\) is the number of pixels per layer. However, this can be optimized by noting that each propagation operator is described by a Toeplitz matrix since the point-spread function in Eq. (3) is invariant with respect to translation in space. Hence, the application of this operator constitutes a convolution operation with the input field. This convolution can be computed by applying the Fourier transform to the input field, followed by multiplication by the Fourier image of the Toeplitz matrix, followed by the inverse Fourier transform. By utilizing the fast Fourier transform algorithm, this operation can be completed in \(\mathcal{O}(N\log N)\) time [42]. Additionally, the transmission operator \(\hat{T}_{l}\) is described by a diagonal matrix and can be evaluated in \(\mathcal{O}(N)\) time. Therefore, calculating the evolution of \(B\) different input fields through a DONN with \(L\) layers and \(N\) pixels per layer has a computational complexity of \(\mathcal{O}(BLN\log N)\), and the backward pass has the same complexity.

### Arbitrary spatial coherence illumination

Using DONNs for real-world applications requires the ability to process incoherent and partially coherent light. We assume quasimonochromatic illumination conditions, which is a good approximation for many cases. These conditions require that the input light is narrowband and its coherence length is much greater than the maximum path length difference between diffractive layers [43]. At the same time, we assume the coherence time to be much shorter than the inverse detection bandwidth, so the detection averages in time over the non-stationary interference pattern.

The spatial coherence of the optical field in the \(l\)-th layer is characterized by the mutual intensity function, which determines the time-averaged correlation of the field at two separate pixels [40; 43]:

\[J_{l}(m,m^{\prime})=\lim_{T\rightarrow\infty}\frac{1}{T}\int_{-T/2}^{T/2}\psi _{l}(m;t)\,\psi_{l}^{*}(m^{\prime};t)\,\mathrm{d}t,\]

where \(T\) is the detection time. This matrix represents an operator

\[\hat{J}_{l}=\lim_{T\rightarrow\infty}\frac{1}{T}\int_{-T/2}^{T/2}|\psi_{l}(t) \rangle\!\langle\psi_{l}(t)|\ \mathrm{d}t, \tag{8}\]

where \(J_{l}(m,m^{\prime})=\langle m|\hat{J}_{l}|m^{\prime}\rangle\). The time-averaged intensity of the field is given by the diagonal elements of the mutual intensity operator, such that

\[I_{l}(m)=J_{l}(m,m). \tag{9}\]

Similar to the evolution of coherent fields through DONNs, the evolution of the mutual intensity operator can be expressed using the transmission and propagation operators. The input mutual intensity operator \(\hat{J}_{\text{in}}\) describes the spatial coherence of the initial field that encodes the input image to be processed by the DONN. The output mutual intensity operator is given by

\[\hat{J}_{\text{out}}=\hat{U}\,\hat{J}_{\text{in}}\,\hat{U}^{\dagger}, \tag{10}\]
where \(\hat{U}\) is the evolution operator of the DONN defined in Eq. (5). The output of the DONN corresponds to the total time-averaged intensity, defined in Eq. (9), incident on the spatial detection regions along the output plane:

\[o(c)=\sum_{n\in\mathcal{D}_{\varepsilon}}J_{\text{out}}(n,n). \tag{11}\]

The evolution of the mutual intensity operator and the corresponding DONN output can be simulated on a computer using Eqs. (10) and (11). Analogous to the previously discussed method, the fast Fourier transform can be leveraged to evaluate the propagation operator \(\hat{P}_{l}\) applied to an arbitrary mutual intensity operator described by an \(N\times N\) matrix, which scales as \(\mathcal{O}(N^{2}\log N)\). The transmission operator \(\hat{T}_{l}\) can similarly be evaluated in \(\mathcal{O}(N^{2})\) time. Hence, simulating the evolution of \(B\) different input fields with arbitrary spatial coherence through a DONN with \(L\) layers and \(N\) pixels per layer has a computational complexity of \(\mathcal{O}(BLN^{2}\log N)\). The backward pass executed during training has the same computational complexity.

### Incoherent illumination

DONNs with fully incoherent input can be treated using the computational approach discussed in the previous section. However, computational costs can be amortised for multiple input examples using the impulse response that characterizes the system.

Fully incoherent input is described by the diagonal mutual intensity operator

\[J_{\text{in}}(m,m^{\prime})=I_{\text{in}}(m)\,\delta_{m,m^{\prime}}, \tag{12}\]

where the time-averaged intensity \(I_{\text{in}}(m)\) encodes information from the \(m\)-th pixel of the input image. We can express the corresponding time-averaged intensity along the output plane, using Eqs. (9) and (10), as

\[I_{\text{out}}(n)=\sum_{m}I_{\text{in}}(m)\cdot\Big{|}\left\langle n|\hat{U}| m\right\rangle\Big{|}^{2}. \tag{13}\]

Here, \(|\left\langle n|\hat{U}|m\right\rangle|^{2}\) corresponds to the intensity at the \(n\)-th pixel in the output plane from point-source illumination at the \(m\)-th input pixel (i.e., the intensity impulse response of the system).

The intensity impulse response of a DONN with \(L\) layers and \(N\) pixels per layer can be determined by calculating the coherent evolution of all \(N\) input pixels through the system, which has a computational complexity of \(\mathcal{O}(LN^{2}\log N)\). This is the same complexity as the previous method using the evolution of the mutual intensity operator. However, once the intensity impulse response is known, the output intensity distributions for \(B\) different incoherent input fields can be calculated using Eq. (13) with a reduced computational complexity of \(\mathcal{O}(BN^{2})\). Thus, the total computational cost of simulating DONN inference with fully incoherent input for multiple input examples can be reduced. This technique can be implemented during training to improve the simulation runtime by using mini-batches of input examples, as calculating the unit responses of the DONN is only required once for each mini-batch. A summary of the computational complexities of simulating DONNs under coherent, arbitrary coherence, and incoherent illumination is shown in Table 1.

### Expressivity of DONNs

The expressive power of DONNs is dependent on the spatial coherence of the input light. Under coherent illumination, DONNs have been shown to outperform linear models [24; 27]. Since the coherent field evolves linearly through the system, this improvement in performance results from the nonlinear intensity measurement of the complex-valued field at the output plane (6), followed by the linear summation of the intensities over the detection regions (7). Thus, DONNs using coherent illumination can be understood as standard neural networks that consist of a complex-valued linear layer with a nonlinear activation function, followed by a real-valued linear layer, as shown in Fig. 1b.

In contrast, under incoherent illumination, the time-averaged output intensity is the sum of the intensity patterns from individual pixel sources in the input plane. Therefore, DONNs with incoherent input illumination cannot perform better than a linear model, as the time-averaged input and output intensity distributions are linearly related according to the intensity impulse response, as shown in Eq. (13).

The improved expressive power of a coherently illuminated DONN arises from the off-diagonal elements in the input mutual intensity operator, which are absent for incoherent light. These elements represent the spatial coherence between two different pixels in the input image. For an arbitrary input mutual intensity operator \(J_{\text{in}}\), the output intensity can be expressed, using Eqs. (9) and (10), as

\[I_{\text{out}}(n)=\sum_{m,m^{\prime}}J_{\text{in}}(m,m^{\prime})\cdot\langle n |\hat{U}|m\rangle\cdot\langle m^{\prime}|\hat{U}^{\dagger}|n\rangle\,. \tag{14}\]

This summation includes off-diagonal elements of the mutual intensity operator, which depend nonlinearly on the input field. Due to this nonlinearity, the performance of DONNs under partially coherent illumination can surpass that with incoherent light, as demonstrated in the following section.

### Performance on MNIST dataset

Using the formalism introduced in the previous section, we trained simulated DONN models to identify handwrit

\begin{table}
\begin{tabular}{l l} \hline \hline Spatial Coherence & Computational Complexity \\ \hline Coherent & \(\mathcal{O}(BLN\log N)\) \\ Arbitrary Coherence & \(\mathcal{O}(BLN^{2}\log N)\) \\ Incoherent & \(\mathcal{O}(BN^{2})+\mathcal{O}(LN^{2}\log N)\) \\ \hline \hline \end{tabular}
\end{table}
Table 1: Computational complexity of modeling the evolution of \(B\) input examples through a DONN with \(L\) layers and \(N\) pixels per layer under different illumination conditions.
ten digits from zero to nine using incoherent, partially coherent, and coherent illumination. The models were trained over 50 epochs using 55,000 images (plus 5,000 images for validation) from the MNIST dataset, each consisting of \(28\times 28\) pixels [44]. The DONNs are composed of five diffractive layers, each with \(100\times 100\) pixels, which modulate the phase of the incident light and are spaced 5 cm apart. Each model was trained using a uniform, normalized optical field incident on the input image of the handwritten digit with a central wavelength of 700 nm. The cross-entropy loss function was used during training to calculate the output error of the model. Each pixel in the diffractive layers has a surface area of \(10\times 10\) um\({}^{2}\), while each pixel in the input pattern is \(30\times 30\) um\({}^{2}\). Each detection region in the output plane, which is associated with a unique digit, has an area of \(250\times 250\) um\({}^{2}\). DONNs trained under coherent and incoherent illumination are illustrated in Fig. 2.

The spatial coherence of the input illumination is given by \(J_{\mathrm{in}}(m,m^{\prime})=\sqrt{I_{\mathrm{in}}(m)I_{\mathrm{in}}^{*}(m^{ \prime})}\,\mu^{\prime_{\mathrm{norm}}(m,m^{\prime})}\), where \(I_{\mathrm{in}}(m)\) is the time-averaged intensity at the \(m\)-th pixel in the input image, \(r_{\mathrm{norm}}(m,m^{\prime})\) is the Euclidean distance between the \(m\)-th and \(m^{\prime}\)-th input pixels, normalized by the pixel pitch, and \(\mu\) quantifies the degree of spatial coherence: \(\mu=1\) for fully coherent and \(\mu=0\) for fully incoherent light. We first trained two DONNs models to process handwritten digits using fully coherent and incoherent illumination. During the training phase, we saved the model parameters that yielded the highest validation accuracy. We then evaluated the performance of these models using

Figure 2: Evolution of input examples in a DONN with five layers trained using coherent and incoherent illumination. **a** The top row shows the phase-only modulation profiles of the five diffractive layers, which were trained to process the MNIST dataset using coherent illumination. In the middle row, the digit zero is fed into the system in the input layer and the time-averaged intensity at each layer is shown. The detection regions at the output layer are indicated using blue and red boxes, where the red box corresponds to the target region. The output of the DONN, which is the total intensity in each detection region, is shown on the right. In the bottom row, the evolution of the digit three is shown. The same scaling for the intensity values is used across each row. **b** The equivalent visualization of a DONN trained using incoherent illumination.
a test set of 10,000 images that were not shown during training. The models trained using coherent and incoherent light achieved test accuracies of 97.54% and 91.17%, respectively. The validation accuracy attained during training, as well as the performance of the models on the test set, are shown in Fig. 3.

We then trained DONNs using input illumination with partially coherent light, where each model was trained to process light with a different degree of spatial coherence. The performance of the models was evaluated, and the test accuracies are shown in Fig. 4a. We also evaluated the robustness of the models by testing their performance using input light with degrees of spatial coherence different from that used during training (Fig. 4b). Not surprisingly, the best performance is achieved when the model is evaluated under the same coherence conditions used during training. However, models trained using incoherent illumination are more robust against changes in the spatial coherence.

## III Discussion

We have demonstrated that the performance of DONNs is dependent on the spatial coherence of the incident illumination. Models using incoherent illumination cannot outperform linear models for information processing tasks. However, as demonstrated in Fig. 4a, the degree of spatial coherence required to achieve optimal performance need not be high: \(\mu\sim 0.6\) means that the mutual coherence between points separated by four pixels is reduced by a factor of 0.13. That is, performance almost at the fully coherent level can be reached even when the transverse coherence length is much less than the size of a whole MNIST digit. This implies that neighboring pixels contain more relevant information for pattern recognition compared to distant pixels. As a result, the DONN model can capture relevant nonlinear relationships in the input data without requiring full spatial coherence between all pixels. In addition to the

Figure 4: Performance of DONN models on the MNIST dataset with varied coherence. **a** Test accuracy achieved with input light of a specified degree of spatial coherence \(\mu\). **b** Accuracies attained by models trained with illumination of spatial coherence \(\mu_{\mathrm{train}}\) but tested using input illumination with spatial coherence \(\mu_{\mathrm{test}}\). The test accuracies along the diagonal correspond to using the same illumination conditions for training and testing, as shown in (a).

Figure 3: DONN training results on the MNIST dataset under coherent and incoherent illumination. **a** Validation accuracy attained by the models at each epoch during training. **b, c** Confusion matrix of the model trained using coherent (b) and incoherent (c) illumination evaluated on the test set.
coherence of the incident light, several other factors affect DONN performance including the training hyperparameters, the axial distances between layers, and the size of the detection regions.

We emphasize that the above relation between the input coherence and DONN expressivity assumes that no further processing of the DONN data is implemented. If, for example, the DONN is followed by an electronic neural network with nonlinear activation layers, the DONN can surpass a linear model even if illuminated incoherently. For example, Rahman _et al._ trained a DONN to classify MNIST by associating two detection regions with each digit and then applying a rational function to compute the network prediction from the intensities of these regions. In this way, the accuracy reached was above that of a linear classifier [39].

Incoherently illuminated DONNs are more broadly applicable to real-world environments, as coherent illumination requires a laser source. However, some degree of coherence can also be achieved by illuminating the object with a distant incoherent source of narrow spatial extent according to the van Cittert -- Zernike theorem [40]. As discussed above, illumination with even a short transverse coherence length can significantly enhance the DONN performance.

The robustness of the system can be improved by training the DONN using various illumination conditions, which could be useful for applications that require DONN operation in different environments. Moreover, the system can be further generalized to operate under a continuum of central frequencies, which has been recently experimentally demonstrated using coherent light [45].

###### Acknowledgements.

 This work is supported by Innovate UK Smart Grant 10043476. MJF is funded by the European Union's Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie grant agreement No. 956071.
