Direct tomography of high-dimensional density matrices for general quantum states of photons
Yiyu Zhou, Jiapeng Zhao, Darrick Hay, Kendrick McGonagle, Robert W. Boyd, Zhimin Shi
DDirect tomography of high-dimensional density matrices for general quantum states
Yiyu Zhou, Jiapeng Zhao, Darrick Hay, Kendrick McGonagle, Robert W. Boyd,
1, 3 and Zhimin Shi ∗ The Institute of Optics, University of Rochester, Rochester, New York 14627, USA Department of Physics, University of South Florida, Tampa, Florida 33620, USA Department of Physics, University of Ottawa, Ottawa, Ontario K1N 6N5, Canada (Dated: February 3, 2021)Quantum state tomography is the conventional method used to characterize density matrices for general quan-tum states. However, the data acquisition time generally scales linearly with the dimension of the Hilbert space,hindering the possibility of dynamic monitoring of a high-dimensional quantum system. Here, we demonstrate adirect tomography protocol to measure density matrices in the spatial domain through the use of a polarization-resolving camera, where the dimension of density matrices can be as large as 580 ×
580 in our experiment. Theuse of the polarization-resolving camera enables parallel measurements in the position and polarization basisand as a result, the data acquisition time of our protocol does not increase with the dimension of the Hilbertspace and is solely determined by the camera exposure time (on the order of 10 milliseconds). Our methodis potentially useful for the real-time monitoring of the dynamics of quantum states and paves the way for thedevelopment of high-dimensional, time-efficient quantum metrology techniques.
Introduction. —The ability to characterize a quantum stateis crucial in quantum technologies, both because it ensuresthat the desired quantum state has been generated and it canbe used to determine the quantum state after interacting with asystem. Quantum state tomography is an established approachto reconstruct a general quantum state (either pure or mixed)through a series of projective measurements performed onidentically prepared states [1–13]. Recently, the concept ofdirect measurement [14] has been established, which can di-rectly be used to read out the complex-valued amplitudes of apure quantum state through a proper sequence of weak [15–25] and strong measurements. The elimination of the com-plicated post-processing procedure of state reconstruction isone of the main advantages of direct measurement methods,allowing it to serve as an alternative metrology technique thatmay greatly reduce experimental complexity.The concept of direct measurement is quickly being ex-tended to the characterization of various quantum systems[26–32]. Nonetheless, one remaining challenge in quantum-state metrology is the limited characterization speed and ef-ficiency for high-dimensional quantum states. Most demon-strated techniques, including direct measurement methods, in-volve either a slow scanning process or a complicated post-processing procedure, where the characterization time scalesunfavorably with the dimension of the quantum system. As aresult, almost all quantum metrology demonstrations to datehave been carried out under stable laboratory conditions, andthe measurement of a high-dimensional quantum state cantake as long as several hours. Compressive sensing has beenimplemented for the tomography of an N -dimensional purestate in the spatial domain, and an experiment with N =19 , has been demonstrated by using only ≈ . N mea-surements [30]. Direct measurement of the density matrix inthe high-dimensional orbital-angular-momentum (OAM) ba-sis has also been reported [31, 32]. However, these meth-ods use single-pixel detectors for data collection and requireperforming a series of measurements via scanning for the re-construction of high-dimensional quantum states. In general, since the number of measurements scales linearly with thedimension of the Hilbert space, the data acquisition time in-evitably increases for high-dimensional quantum states, hin-dering the possibility of real-time monitoring of dynamicquantum systems. While the recently proposed auxiliaryHilbert space tomography [33] can reduce the measurementcomplexity for density matrix characterization, this methodis only applicable to OAM states and thus exhibits a lim-ited range of application. In the following, we introduce ascan-free direct tomography protocol that can measure thecomplex-valued high-dimensional density matrix for mixedphoton states by using a polarization-resolving camera. Thedata acquisition time does not increase with the dimension ofthe Hilbert space, and the maximum dimension allowed byour protocol is only limited by the pixel count of the detectorarray. Direct tomography protocol. —The density matrix can berepresented as an incoherent mixture of pure states, which canbe expressed as ˆ ρ = (cid:88) k p k | ψ k (cid:105) (cid:104) ψ k | , (1)where | ψ k (cid:105) is the pure quantum state normalized as (cid:104) ψ k | ψ k (cid:105) = 1 , and p k is the probability coefficient normalizedas (cid:80) k p k = 1 . The element of density matrix in the positionbasis can be computed as ρ ( x , y , x , y ) = (cid:104) x , y | ˆ ρ | x , y (cid:105) = (cid:88) k p k (cid:104) x , y | ψ k (cid:105) (cid:104) ψ k | x , y (cid:105) = (cid:88) k p k ψ k ( x , y ) ψ ∗ k ( x , y ) , (2)where | x, y (cid:105) denotes the position eigenstate located at ( x, y ) .In our experiment, we assume that the traverse profile ofthe quantum state has only a one-dimensional (1D) vari-ation along the x axis and is invariant along the y axis.Therefore, we have ψ k ( x , y ) = ψ k ( x ) is independent a r X i v : . [ qu a n t - ph ] F e b of y , and the density matrix element can be simplified as ρ ( x , y , x , y ) = ρ ( x , x ) , and ρ ( x , x ) = (cid:88) k p k ψ k ( x ) ψ ∗ k ( x ) , (3)which is the quantity to be measured. It is worth notingthat Eq. (3) is reminiscent of the mutual coherence functionin classical optics [34, 35]. We use the polarization as thethe pointer state [26] which is prepared in the diagonal po-larization state | D (cid:105) = ( | H (cid:105) + | V (cid:105) ) / √ , where | H (cid:105) and | V (cid:105) denote the horizontal and vertical polarization state, respec-tively. Therefore, the full initial density matrix can be writtenas ˆ ρ i = ˆ ρ ⊗ | D (cid:105) (cid:104) D | . Our direct tomography protocol entailsintroducing a 90° beam rotation for the horizontally polarizedbeam while leaving the vertically polarized beam unchanged.This polarization-sensitive beam rotation can be described bya unitary transformation as [32] ˆU = ˆT( π/ ⊗ | H (cid:105) (cid:104) H | + ˆT(0) ⊗ | V (cid:105) (cid:104) V | , (4)where ˆT( θ ) = exp (cid:16) iθ ˆ (cid:96) (cid:17) is the rotation operator, ˆ (cid:96) is the or-bital angular momentum operator about the optical axis, andthe effect of rotation operator on the position eigenstate can bewritten as ˆT( θ ) | x, y (cid:105) = | x cos θ − y sin θ, x sin θ + y cos θ (cid:105) .The final density matrix after this unitary transformation canbe represented as ˆ ρ f = ˆU † ˆ ρ i ˆU . The projective measurements[36, 37] we propose to perform can be represented by the fol-lowing projectors: ˆ π D = | x, y (cid:105) (cid:104) x, y | ⊗ | D (cid:105) (cid:104) D | , ˆ π A = | x, y (cid:105) (cid:104) x, y | ⊗ | A (cid:105) (cid:104) A | , ˆ π R = | x, y (cid:105) (cid:104) x, y | ⊗ | R (cid:105) (cid:104) R | , ˆ π L = | x, y (cid:105) (cid:104) x, y | ⊗ | L (cid:105) (cid:104) L | , (5)where | A (cid:105) = ( | H (cid:105) − | V (cid:105) ) / √ is the anti-diagonal polariza-tion state, | L (cid:105) = ( | H (cid:105) + i | V (cid:105) ) / √ is the left-handed circularpolarization state, and | R (cid:105) = ( | H (cid:105) − i | V (cid:105) ) / √ is the right-handed circular polarization state. Therefore, the expectationvalue of these projectors are found to be Γ D ( x, y ) = Tr[ ˆ π D ˆU † ˆ ρ i ˆU]= 14 ( ρ ( − y, − y ) + ρ ( x, x ) + 2Re[ ρ ( − y, x )]) , Γ A ( x, y ) = Tr[ ˆ π A ˆU † ˆ ρ i ˆU]= 14 ( ρ ( − y, − y ) + ρ ( x, x ) − ρ ( − y, x )]) , Γ R ( x, y ) = Tr[ ˆ π R ˆU † ˆ ρ i ˆU]= 14 ( ρ ( − y, − y ) + ρ ( x, x ) + 2Im[ ρ ( − y, x )]) , Γ L ( x, y ) = Tr[ ˆ π L ˆU † ˆ ρ i ˆU]= 14 ( ρ ( − y, − y ) + ρ ( x, x ) − ρ ( − y, x )]) . (6)Using the above equations, the density matrix can be experi- SLM LensIris PolarizerHWP
P B S D o v e p r i s m QWP PolarCamPBS H e N e l a s e r xy FIG. 1. The schematic of the experimental setup. A Dove prismis used to rotate the horizontally polarized beam by 90°. A 1DHermite-Gauss state HG ( x ) is used as an example to visualize thebeam rotation. SLM: spatial light modulator. HWP: half-wave plate.PBS: polarizing beamsplitter. QWP: quarter-wave plate. PolarCam:polarization-resolving camera. mentally reconstructed as ρ exp0 ( x , x ) = Γ D ( x , − x ) − Γ A ( x , − x )+ i (Γ R ( x , − x ) − Γ L ( x , − x )) . (7)It can be seen that the density matrix can be directly recon-structed without using any complicated algorithm. In additionto the reconstruction of the density matrix ˆ ρ exp0 , it is also de-sirable to be able to reconstruct the pure states | ψ k (cid:105) . In orderto reconstruct the pure states, we use singular value decompo-sition [38]. The reconstruction can be unique if the pure statesare mutually orthogonal. For a square and Hermitian densitymatrix ˆ ρ exp0 , it can always be decomposed as [39] ˆ ρ exp0 = ˆMˆS ˆM † , (8)where ˆM is a unitary matrix and ˆS is a real-valued diagonalmatrix whose diagonal elements S kk are the singular valuesof ˆ ρ exp0 . It can be readily seen that (cid:104) x | ˆ ρ exp0 | x (cid:105) = (cid:88) k S kk (cid:104) x | ˆM | k (cid:105) (cid:104) k | ˆM † | x (cid:105) . (9)Comparing Eq. (9) with Eq. (2), one can find that p exp k = S kk , (cid:104) x | ψ exp k (cid:105) = (cid:104) x | ˆM | k (cid:105) . (10)As one can see, singular value decomposition can be used as atool to decompose a density matrix into an incoherent mixtureof pure states, which can be efficiently implemented by estab-lished numerical algorithms [40]. It is worth noting that thesingular value decomposition discussed here is reminiscent ofthe coherent mode decomposition in optical coherence theory[41, 42]. Γ D ( x ,- x ) Γ A ( x ,- x ) Γ R ( x ,- x ) Γ L ( x ,- x ) N o r m a li z ed i n t en s i t y Experiment Re(ρ )1 5801580 Theory Re(ρ )1 5801580 Experiment Im(ρ )1 5801580 Theory Im(ρ )1 5801580 -3 -2×10 -3 A m p li t ude x (pixel) x ( p i x e l ) x (pixel) x (pixel) x (pixel) x ( p i x e l ) x ( p i x e l ) x ( p i x e l ) x (pixel) x (pixel) x (pixel) x (pixel) x ( p i x e l ) x ( p i x e l ) x ( p i x e l ) x ( p i x e l ) ab exp exp FIG. 2. Experimental results for the phase-only states. (a) The images acquired by the PolarCam. (b) The real and imaginary part of thereconstructed density matrix. The trace distance between the theoretical density matrix and the experimentally measured density matrix is14.2% ± Γ D ( x ,- x ) Γ A ( x ,- x ) Γ R ( x ,- x ) Γ L ( x ,- x ) N o r m a li z ed i n t en s i t y Experiment Re(ρ )1 5801580 Theory Re(ρ )1 5801580 Experiment Im(ρ )1 5801580 Theory Im(ρ )1 5801580 -3 -5×10 -3 A m p li t ude x (pixel) x ( p i x e l ) x (pixel) x (pixel) x (pixel) x ( p i x e l ) x ( p i x e l ) x ( p i x e l ) x (pixel) x (pixel) x (pixel) x (pixel) x ( p i x e l ) x ( p i x e l ) x ( p i x e l ) x ( p i x e l ) ab exp exp FIG. 3. Experimental results for the HG states. (a) The images acquired by the PolarCam. (b) The real and imaginary part of the reconstructeddensity matrix. The trace distance between the theoretical density matrix and the experimentally measured density matrix is 19.0% ± Experiment. —The experimental setup to implement the di-rect tomography protocol is shown in Fig. 1. A 633 nm HeNelaser is used as the source of photons. The light beam isspatially filtered and attenuated before it illuminates a spa-tial light modulator (SLM, Pluto 2 VIS-020, Holoeye). A se-ries of computer-generated phase-only holograms [43] is dis-played onto the SLM to generate the quantum states of inter-est. Mixed states can be generated by switching the hologramon the SLM and by incoherently mixing the intensity imagesacquired by the camera [44]. An iris is used to pass the firstdiffraction order of light coming off the SLM while blockingall other diffraction orders. A polarizer and a half-wave plate (HWP) are used to generate the diagonal polarization state | D (cid:105) . To implement the unitary transformation ˆU [cf. Eq. (4)],a polarizing beamsplitter (PBS) is used to separate the hori-zontally and vertically polarized beam. A Dove prism is ap-plied to geometrically rotate the horizontally polarized beamby 90°. A second PBS is used to recombine the two beamsand thus completes the implementation of ˆU . A 45°-orientedquarter-wave plate (QWP) and a polarization-resolving cam-era (PolarCam, BFS-U3-51S5P-C, FLIR) are used to performall the required projective measurements in a single shot. ThePolarCam has micro-sized polarizers (oriented to 0°, 45°, 90°,and 135°) deposited on the camera sensors and thus allows A m p li t ude A m p li t ude Re( p ψ ( x ) ) Im( p ψ ( x ) ) x (pixel) x (pixel) b A m p li t ude A m p li t ude Re( p ψ ( x ) ) Im( p ψ ( x ) ) x (pixel) x (pixel) a Re( p ψ ( x ) ) Im( p ψ ( x ) ) A m p li t ude x (pixel) A m p li t ude x (pixel) c State fidelity |⟨ ψ | ψ ⟩| =92.6%±2.1% exptheory State fidelity |⟨ ψ | ψ ⟩| =93.1%±1.7% exptheory State fidelity |⟨ ψ | ψ ⟩| =96.5%±1.0% exptheory FIG. 4. The reconstructed phase-only quantum state for (a) | ψ (cid:105) , (b) | ψ (cid:105) , and (c) | ψ (cid:105) , respectively. The real (imaginary) part is shownin the left (right) panel. The standard deviation of the experimentaldata is denoted by the line width, which is generally too small to bevisible. The fidelity of each reconstructed state is shown at the top ofeach corresponding subfigure. for the detection of four different polarization states simulta-neously. The camera exposure time is approximately 10 mil-liseconds depending on the intensity of the generated states.The QWP and the PolarCam jointly enable the projective mea-surements proposed in Eq. (5). The image on the camera hasa size of × pixels, and thus the dimensionality of thequantum states in our experiment is N = 580 . The pixel sizeof the camera is 3.45 µ m.In our experiment, we prepare a mixed state consisting ofthree mutually orthogonal pure states | ψ k (cid:105) with k = 1 , , .More specifically, as our first demonstration, the pure statesused to construct the density matrix are p = 0 . , (cid:104) x | ψ (cid:105) = e i . πx/a ,p = 0 . , (cid:104) x | ψ (cid:105) = e iπ [ − . x/a ) +4 . x/a )] ,p = 0 . , (cid:104) x | ψ (cid:105) = e iπ [ − . x/a ) − x/a ] , (11)where − a/ < x ≤ a/ is the discretized position, and a =2 mm is the size of the beam. These states are referred to thephase-only quantum states henceforth.As another test of our protocol, we use the 1D Hermite- State fidelity |⟨ ψ | ψ ⟩| =92.4%±0.7% exptheory A m p li t ude A m p li t ude Re( p ψ ( x ) ) Im( p ψ ( x ) ) bac A m p li t ude A m p li t ude Re( p ψ ( x ) ) Im( p ψ ( x ) ) x (pixel) x (pixel) Re( p ψ ( x ) ) Im( p ψ ( x ) ) A m p li t ude x (pixel) A m p li t ude x (pixel) ExperimentTheory x (pixel) x (pixel)State fidelity |⟨ ψ | ψ ⟩| =92.9%±0.3% exptheory State fidelity |⟨ ψ | ψ ⟩| =96.8%±0.7% exptheory FIG. 5. The reconstructed HG quantum state for (a) | ψ (cid:105) , (b) | ψ (cid:105) ,and (c) | ψ (cid:105) , respectively. The real (imaginary) part is shown in theleft (right) panel. The standard deviation of the experimental data isdenoted by the line width, which is generally too small to be visible.The fidelity of each state reconstructed is shown at the top of eachcorresponding subfigure. Gauss (HG) states to construct the mixed state: HG m ( x ) = (cid:18) πw (cid:19) √ m m ! × H m (cid:32) √ xw (cid:33) exp (cid:18) − x w (cid:19) , (12)where H m ( · ) is the Hermite polynomial of order m [45], and w = 0 . a is the beam waist radius. The HG states used inthe experiment are p = 0 . , (cid:104) x | ψ (cid:105) = HG ( x ) ,p = 0 . , (cid:104) x | ψ (cid:105) = HG ( x ) ,p = 0 . , (cid:104) x | ψ (cid:105) = HG ( x ) . (13)The images acquired by the PolarCam for the phase-onlystates are shown in Fig. 2(a). We apply a digital low-passGaussian spatial filter to process these images in order to re-move the undesirable fringes caused by dusts and glass filminterference. The density matrix can be directly reconstructedbased on these data by using Eq. (7). Due to the experimen-tal errors (e.g., misalignments, noises, imperfect mode gener-ation fidelity, etc.), the experimentally reconstructed densitymatrix ˆ ρ exp0 may not be strictly Hermitian. Hence, we imple-ment ˆ ρ exp0 → (ˆ ρ exp0 + ˆ ρ exp † ) / to guarantee the Hermiticity ofthe density matrix, and the results are shown in Fig. 2(b). Toquantify the accuracy of our protocol, we calculate the tracedistance between the ideal density matrix ˆ ρ and the experi-mentally measured density matrix ˆ ρ exp0 as follows [38]:Trace distance = 12 (cid:12)(cid:12)(cid:12)(cid:12) Tr[ (cid:113) ( ˆ ρ − ˆ ρ exp0 )( ˆ ρ − ˆ ρ exp0 ) † ] (cid:12)(cid:12)(cid:12)(cid:12) , (14)and the trace distance for the phase-only states is calculatedto be 14.2% ± ± | ψ exp k (cid:105) ,we compute its fidelity as (cid:12)(cid:12)(cid:12)(cid:68) ψ theory k (cid:12)(cid:12)(cid:12) ψ exp k (cid:69)(cid:12)(cid:12)(cid:12) . The fidelityfor each reconstructed state is shown at the top of each cor-responding subfigure. It can be seen that the fidelity of stateis always higher than 90%. In our experiment, we attributethe nonzero trace distance primarily to the imperfect spatialmode generation and the misalignment of the polarization-sensitive beam rotator. As a consequence, the reconstructeddensity matrix might be unphysical due to the possible lackof Hermiticity and positive semi-definiteness [46]. However,we notice that the standard maximum-likelihood-estimation-based routine for recovering a physical density matrix [46] isnot readily applicable to our experiment, because it requiresthe minimization of a likelihood function with N = 336 , independent parameters. This task can potentially be accom-plished by using machine learning algorithms [47] and is sub-ject to future study. In our experiment, we assume the trans-verse profile of the field has only a 1D variation along the x axis [see Eq. (3)]. Although our protocol cannot be directlyapplied to a general two-dimensional (2D) spatial field, it ispossible to use a SLM to realize the 2D-to-1D beam reshap-ing (similar to the conversion from a 2D matrix to a 1D vector)and then apply our protocol subsequently.It is worthwhile to mention that many quantum-state-metrology techniques available today use single-pixel detec-tors, such as single-photon avalanche diodes (SPADs). How-ever, advances in detector development in the past few yearshave led to many options for the use of high-performancedetector arrays. It has been demonstrated experimentallythat SPAD arrays [48, 49], commercial cooled CCD cam-eras [10, 23], electron-multiplying CCD cameras [50, 51] andintensified CCD cameras [52–57] can all operate at single-photon levels. All pixels in the array can detect incomingphotons simultaneously during the exposure time. Compar-ing with raster scanning techniques using a single-pixel de- tector, the parallel measurement via a N -pixel detector ar-ray is generally N times faster. The emergence of these newsingle-photon-sensitive detector arrays can readily enable usto develop new and more efficient quantum state metrologyprotocols that can best utilize the availability of parallel mea-surements at single photon levels. Conclusion. —In this work, we demonstrated a direct to-mography protocol that can efficiently characterize a high-dimensional density matrix in the spatial domain for generalquantum states, and the data acquisition time is independentof the dimension of the Hilbert space. Polarization was usedas the pointer state, and a polarization-sensitive beam rotationwas applied to enable the detection of off-diagonal elements inthe density matrix. Two different mixed states were preparedand characterized with a high fidelity in our demonstration.Singular value decomposition was implemented to reconstructthe pure states that constitute the prepared mixed state, whichcan potentially be useful for the analysis of spatially incoher-ent fields. Our protocol can readily scale to higher dimensionswithout increasing the data acquisition time. Furthermore, themaximum measurable dimension is limited only by the cam-era pixel number. We anticipate that our protocol can inspirethe development of high-dimensional, time-efficient quantummetrology techniques and can be used as a powerful tool forthe experimental study of the spatial mutual coherence func-tion of optical fields, which play an important role in super-resolution imaging [58] and the optical coherence theory [34].This work is supported by the U.S. Office of Naval Re-search (N00014-17-1-2443 and N00014-20-1-2558). In addi-tion, R.W.B. acknowledges support from Canada ExcellenceResearch Chairs Program and Natural Sciences and Engineer-ing Research Council of Canada. ∗ [email protected][1] A. G. White, D. F. V. James, W. J. Munro, and P. G. Kwiat,Phys. Rev. A , 012301 (2001).[2] J. Itatani, J. Levesque, D. Zeidler, H. Niikura, H. P´epin, J. C.Kieffer, P. B. Corkum, and D. M. Villeneuve, Nature , 867(2004).[3] K. J. Resch, P. Walther, and A. Zeilinger, Phys. Rev. Lett. ,070402 (2005).[4] J. S¨oderholm, G. Bj¨ork, A. B. Klimov, L. L. S´anchez-Soto, andG. Leuchs, New J. Phys. , 115014 (2012).[5] D. Sych, J. ˇReh´aˇcek, Z. Hradil, G. Leuchs, and L. L. S´anchez-Soto, Phys. Rev. A , 052123 (2012).[6] A. I. Lvovsky, H. Hansen, T. Aichele, O. Benson, J. Mlynek,and S. Schiller, Phys. Rev. Lett. , 050402 (2001).[7] M. Beck, Phys. Rev. Lett. , 5748 (2000).[8] M. Paris and J. Rehacek, Quantum state estimation (SpringerScience & Business Media, 2004).[9] M. Beck, C. Dorrer, and I. A. Walmsley, Phys. Rev. Lett. ,253601 (2001).[10] A. M. Dawes, M. Beck, and K. Banaszek, Phys. Rev. A ,032102 (2003).[11] B. J. Smith, B. Killett, M. G. Raymer, I. A. Walmsley, andK. Banaszek, Opt. Lett. , 3365 (2005). [12] D. T. Smithey, M. Beck, M. G. Raymer, and A. Faridani, Phys.Rev. Lett. , 1244 (1993).[13] A. I. Lvovsky and M. G. Raymer, Rev. Mod. Phys. , 299(2009).[14] J. S. Lundeen, B. Sutherland, A. Patel, C. Stewart, and C. Bam-ber, Nature , 188 (2011).[15] Y. Aharonov, D. Z. Albert, and L. Vaidman, Phys. Rev. Lett. ,1351 (1988).[16] I. M. Duck, P. M. Stevenson, and E. C. G. Sudarshan, Phys.Rev. D , 2112 (1989).[17] N. W. M. Ritchie, J. G. Story, and R. G. Hulet, Phys. Rev. Lett. , 1107 (1991).[18] L. M. Johansen, Phys. Rev. Lett. , 120402 (2004).[19] O. Hosten and P. Kwiat, Science , 787 (2004).[20] D. R. Solli, C. F. McCormick, R. Y. Chiao, S. Popescu, andJ. M. Hickmann, Phys. Rev. Lett. , 043601 (2004).[21] P. B. Dixon, D. J. Starling, A. N. Jordan, and J. C. Howell, Phys.Rev. Lett. , 173601 (2009).[22] A. Feizpour, X. Xing, and A. M. Steinberg, Phys. Rev. Lett. , 133603 (2011).[23] S. Kocsis, B. Braverman, S. Ravets, M. J. Stevens, R. P. Mirin,L. K. Shalm, and A. M. Steinberg, Science , 1170 (2011).[24] J. Dressel, M. Malik, F. M. Miatto, A. N. Jordan, and R. W.Boyd, Rev. Mod. Phys. , 307 (2014).[25] C. Bamber and J. S. Lundeen, Phys. Rev. Lett. , 070405(2014).[26] G. S. Thekkadath, L. Giner, Y. Chalich, M. J. Horton, J. Banker,and J. S. Lundeen, Phys. Rev. Lett. , 120401 (2016).[27] J. S. Lundeen and C. Bamber, Phys. Rev. Lett. , 070402(2012).[28] W. Shengjun, Sci. Rep. , 1193 (2013).[29] J. Z. Salvail, M. Agnew, A. S. Johnson, E. Bolduc, J. Leach,and R. W. Boyd, Nat. Photon. , 316 (2013).[30] M. Mirhosseini, O. S. Maga˜na Loaiza, S. M. Hashemi Rafsan-jani, and R. W. Boyd, Phys. Rev. Lett. , 090402 (2014).[31] M. Malik, M. Mirhosseini, M. P. J. Lavery, J. Leach, M. J. Pad-gett, and R. W. Boyd, Nat. Commun. , 3115 (2014).[32] M. Mirhosseini, O. S. Maga˜na Loaiza, C. Chen, S. M.Hashemi Rafsanjani, and R. W. Boyd, Phys. Rev. Lett. ,130402 (2016).[33] R. Liu, J. Long, P. Zhang, R. E. Lake, H. Gao, D. P. Pappas, andF. Li, arXiv:1908.00577 (2019).[34] L. Mandel and E. Wolf, Optical coherence and quantum optics (Cambridge University Press, Cambridge, 1995).[35] Z. Hradil, J. ˇReh´aˇcek, and L. L. S´anchez-Soto, Phys. Rev. Lett. , 010401 (2010).[36] Z. Shi, M. Mirhosseini, J. Margiewicz, M. Malik, F. Rivera,Z. Zhu, and R. W. Boyd, Optica , 388 (2015).[37] Z. Zhu, D. Hay, Y. Zhou, A. Fyffe, B. Kantor, G. S. Agarwal,R. W. Boyd, and Z. Shi, Phys. Rev. Applied , 034036 (2019).[38] M. A. Nielsen and I. L. Chuang, Quantum Computationand Quantum Information (Cambridge University Press, Cam-bridge, 2010).[39] L. N. Trefethen and D. Bau III,
Numerical linear algebra (SIAM, Philadelphia, 1997).[40] V. Klema and A. Laub, IEEE Trans. Automat. Contr. , 164(1980).[41] E. Wolf, J. Opt. Soc. Am. , 343 (1982).[42] E. Wolf, Opt. Commun. , 3 (1981).[43] V. Arriz´on, U. Ruiz, R. Carrada, and L. A. Gonz´alez, J. Opt.Soc. Am. A , 3500 (2007).[44] B. Rodenburg, M. Mirhosseini, O. S. Maga˜na Loaiza, and R. W.Boyd, J. Opt. Soc. Am. B , A51 (2014).[45] O. Svelto and D. C. Hanna, Principles of lasers (Springer, NewYork, 2010).[46] D. F. V. James, P. G. Kwiat, W. J. Munro, and A. G. White,Phys. Rev. A , 052312 (2001).[47] G. Torlai, G. Mazzola, J. Carrasquilla, M. Troyer, R. Melko,and G. Carleo, Nat. Phys. , 447 (2018).[48] E. Charbon, Philos. Trans. R. Soc. A , 20130100 (2014).[49] G. Gariepy, N. Krstajic, C. L. Robert Henderson, R. R. Thom-son, G. S. Buller, B. Heshmat, R. Raskar, J. Leach, and D. Fac-cio, Nat. Commun. , 6021 (2015).[50] M. Edgar, D. Tasca, F. Izdebski, R. Warburton, J. Leach, M. Ag-new, G. Buller, R. Boyd, and M. Padgett, Nat. Commun. , 984(2012).[51] G. B. Lemos, V. Borish, G. D. Cole, S. Ramelow, R. Lap-kiewicz, and A. Zeilinger, Nature , 409 (2014).[52] D. S. Tasca, R. S. Aspden, P. A. Morris, G. Anderson, R. W.Boyd, and M. J. Padgett, Opt. Express , 30460 (2013).[53] R. S. Aspden, D. S. Tasca, R. W. Boyd, and M. J. Padgett, NewJ. Phys. , 073032 (2013).[54] R. Fickler, M. Krenn, R. Lapkiewicz, S. Ramelow, andA. Zeilinger, Sci. Rep. , 1914 (2013).[55] R. Machulka, O. Haderka, J. Peˇrina, M. Lamperti, A. Allevi,and M. Bondani, Opt. Express , 13374 (2014).[56] P. A. Morris, R. S. Aspden, J. E. Bell, R. W. Boyd, and M. J.Padgett, Nat. Commun. , 5931 (2015).[57] R. Chrapkiewicz, M. Jachura, K. Banaszek, and W. Wasilewski,Nat. Photon. , 576 (2016).[58] W. Larson and B. E. Saleh, Optica5